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ABSTRACT 


Analytic  solutions  to  the  thermal  radiative  transfer  equation  for 
one  and  two-dimensional  cylindrical  geometries  for  an  isothermal  cloud 
cylinder  are  derived  using  a  truncated  spherical  harmonics  expansion 
for  the  intensity  and  scattering  phase  function.  Higher  order  approx¬ 
imations  are  generated  by  utilizing  the  recursion  properties  of  the 
modified  Bessel  functions.  Computational  results  are  presented  and 
physically  discussed  for  a  hypothetical  atmosphere  in  the  infrared 
10  urn  wavelength  and  in  a  realistic  atmosphere  in  the  microwave 
50-60  GHz  frequencies.  Comparisons  are  made  with  a  discrete-ordinate 
plane-parallel  model  and  finite  cubic  and  rectangular  models  in  the 
infrared  region  and  with  a  discrete-ordinate  plane-parallel  model  in 
the  microwave  region.  Results  demonstrate  that  finite  clouds  have 
emissivities  less  than  those  of  plane-parallel  clouds  and  that  opti¬ 
cally  thick  cylindrical  clouds  should  not  be  considered  as  black 
clouds.  Further,  the  effect  of  the  angle  of  observation  of  the  satel¬ 
lite  on  remote  sensing  of  a  cloud  pattern  is  demonstrated  to  be 
significant.  Case  studies  are  presented  utilizing  Defence  Meteorologi¬ 
cal  Satellite  Program  microwave  sounder  data  along  with  standard 
meteorological  data.  Calculations  for  nonprecipitating  and  clear 
cases  were  similar  and  in  good  agreement  with  observations.  ^Precipi¬ 
tating  cloud  cases  were  simulated  well  by  the  cylinder  model, while 
the  plane-parallel  model  failed  in  cases  with  higher  rainfall  rates. 


Finally,  temperature  profile  retrievals  were  performed  using  the  Air 
Force  Global  Weather  Central  statistical  method.  The  method  was  shown 
to  fail  in  case  studies  where  there  was  a  precipitating  cloud  in  the 
field  of  view. 
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CHAPTER  1 


INTRODUCTION 

In  recent  years,  one  of  the  unresolved  problems  in  the  fields  of 
atmospheric  radiation  and  satellite  sensing  has  been  the  transfer  of 
radiation  in  finite  cloud  layers  in  conjunction  with  the  interpreta¬ 
tion  of  satellite  sounding  data.  Additional  problems  include  the 
development  of  sounding  techniques  for  cloud  inference,  and  the  under¬ 
standing  of  the  generation  of  heating  and  cooling  in  the  atmosphere. 
There  are  many  circumstances  in  the  Earth's  atmosphere  by  which 
finite  clouds  and  cloud  bands  are  formed.  Inference  of  the  composi¬ 
tion  and  structure  of  clouds  and  the  atmosphere  utilizing  satellite 
sounding  techniques  requires  a  comprehensive  understanding  of  the 
interaction  of  these  finite  cloud  layers  with  the  radiation  field  of 
atmosphere. 

Several  recent  studies  have  been  carried  out  to  investigate  the 
transfer  of  radiation  in  finite  clouds  using  cartesian  coordinates. 
McKee  and  Cox  (1976),  Davies  (1978),  and  Davis  et.  al.  (1979)  re¬ 
ported  the  significant  effect  of  the  cubic  geometry  on  the  reflection, 
transmission  and  absorption  of  solar  radiation.  Weinman  and  Davies 
(1978)  investigated  the  cubic  geometry  influence  on  the  transfer  of 
microwave  radiation  at  a  frequency  of  37  GHz.  More  recently,  Liou 
and  Ou  (1979),  using  an  infrared  window  wavelength  of  10  pm,  found 
that  finite  geometry  significantly  affects  cloud  emissivity  and  infra¬ 
red  cooling  rates.  None  of  these  studies,  however,  have  applied  the 
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finite  cloud  transfer  model  to  realistic  atmospheres  which  include  the 
absorption  contribution  of  atmospheric  gases  and  actual  satellite 
sounding  problems  in  which  observed  data  are  available.  Moreover, 
transfer  of  radiation  in  coordinate  systems  other  than  rectangular  has 
not  been  theoretically  explored  and  investigated  in  atmospheric  remote 
sounding  and  radiative  transfer  problems. 

The  first  objective  of  this  study  is  to  develop  a  theoretical 
finite  cloud  model  using  a  cylindrical  geometry  applicable  to  fre¬ 
quently  occurring  cloud  and  precipitation  bands  in  the  microwave  and 
infrared  radiation  regions  and  to  apply  the  cylindrical  model  to  real¬ 
istic  atmospheric  conditions.  The  second  objective  is  to  modify  the 
cylindrical  cloud/precipitation  model  to  interpret  the  observed 
brightness  temperatures  sensed  by  the  DMSP  (Defense  Meteorological 
Satellite  Program)  microwave  sounder  (SSM/T)  which  was  launched  in 
June,  1979,  for  carefully  selected  clear,  cloudy  and  precipitating 
cases.  In  connection  with  these  cases  and  depending  on  the  availabil¬ 
ity  of  radiosonde  temperature  profiles,  temperature  retrieval  exer¬ 
cises  based  on  the  observed  brightness  temperatures  are  carried  out 
utilizing  a  statistical  method. 

The  initial  impetus  for  radiative  transfer  in  cylindrical 
geometry  was  provided  by  Mark  (1945)  in  the  solution  in  the  neutron 
transport  or  Boltzman  equation.  From  Mark’s  work,  investigation 
was  begun  in  the  neutron  transport  field  in  a  number  of  different 
geometries.  The  spherical  harmonics  approach  for  the  Boltzman 
equation  in  a  one-dimensional  cylindrical  medium  was  investigated  by 
Davison  (1957)  and  Kofink  (1959).  Heasley  and  Warming  (1966)  applied 
a  moment  method  to  the  one-dimensional  cylindrical  problem  in 
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radiative  transfer,  and  Heasley  (1977)  formulated  the  transfer  equa¬ 
tion  as  a  second-order  differential  equation  resulting  in  a  set  of 
tridiagonal  difference  equations.  The  work  presented  in  this  disser¬ 
tation  more  directly  follows  the  lead  presented  by  Kofink  (1959)  in  the 
neutron  transport  field,  applying  a  modification  of  his  method  to  one¬ 
dimensional  cylindrical  radiative  transfer.  Radiative  transfer  in  a 
two-dimensional  cylindrical  medium  is  then  developed  for  various  orders 
of  approximation. 

In  Chapter  2,  the  detailed  characteristics  of  the  DMSP  microwave 
radiometer  (SSM/T)  are  presented.  Chapter  3  is  devoted  to  an  analysis 
of  the  optical  parameters  peculiar  to  the  frequencies  used  by  the 
SSM/T.  In  Chapter  4,  thermal  radiative  transfer  in  plane-parallel 
atmospheres  utilizing  the  discrete-ordinate  method  is  reviewed  briefly. 
Chapter  5  reports  the  theoretical  development  for  radiative  transfer 
in  cylindrical  media,  utilizing  a  spherical  harmonics  approach.  Chap¬ 
ter  6  presents  computational  results  utilizing  the  two-dimensional 
cylindrical  model.  Comparisons  are  made  with  results  presented  by 
Liou  and  Ou  (1979)  using  the  infrared  10  urn  window  region  for  cubic 
and  rectangular  cloud  models.  Comparisons  are  also  made  with  results 
reported  by  Nipko  (1979)  using  the  microwave  60  GHz  SSM/T  frequencies 
for  the  plane-parallel  model.  In  Chapter  7,  data  sets,  which  include 
the  SSM/T  microwave  sounder  data  and  conventional  meteorological  data, 
are  chosen  for  two  days.  By  employing  the  temperature  and  humidity 
profiles  and/or  cloud  information,  radiative  transfer  calculations  are 
carried  out  and  compared  with  observed  satellite  data.  Finally,  tem¬ 
perature  profile  retrievals  are  performed  using  a  statistical  method. 


CHAPTER  2 


DMSP  BLOCK  50  SSM/T  SOUNDER  CHARACTERISTICS 

The  SSM/T  microwave  scanning  radiometer  was  launched  in  June, 

1979,  in  a  sun -synchronous  polar  orbit,  by  the  U.  S.  Ai  r  Force  as 
part  of  the  Defense  Meteorological  Satellite  Program  (DMSP)  Block  5D 
package.  Other  equipment  of  meteorological  interest  on  board  the 
spacecraft  include  visual  and  infrared  imagery  channels  and  a  scanning 
infrared  spectroradiometer.  The  microwave  instrument  scans  across  the 
track  of  the  satellite  (Figure  1),  sensing  a  near  circle  of  174  km 
diameter  at  nadir,  and  an  ellipse  with  a  major  axis  of  304  km  and 
minor  axis  of  approximately  213  km  at  maximum  scan  angle  (Figure  2) 
of  36°  from  nadir.  The  scan  parameter  information  in  Table  1  may  be 
used  in  conjunction  with  Figure  2  to  visualize  how  the  sounder  "sees" 
the  earth-atmosphere  system  as  it  moves  along  its  track.  Seven  opera¬ 
tional  frequencies  (Table  2)  were  chosen  in  the  vicinity  of  a  strong 
oxygen  absorption  band  near  60  GHz,  so  that  the  weighting  functions 
associated  with  the  various  channels  peak  at  different  altitudes  with¬ 
in  the  Earth's  atmosphere,  thus  providing  sufficient  vertical  resolu¬ 
tion  for  temperature  profile  retrieval.  Channel  1,  50.50  GHz,  is  a 
“window"  channel  and  senses  at  or  near  the  Earth's  surface.  Channels 
2,  3  and  4  (53.20,  54.35  and  54.90  GHz)  have  weighting  functions  that 
peak  at  increasing  heights  in  the  troposphere.  The  weighting  functions 
of  Channels  5,  6  and  7  peak  in  the  stratosphere  with  Channel  5  peaking 
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Figure  I.  SSM/T  subtrack  scan  patten.. 
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Hot  Load 
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Figure  2.  SSM/T  Scan  Geometry.  Source:  Rigone  and  Stogryn  (  1977) 


Table  1.  Scan  Parameters 


Scan  Type . 

Cross-Track  Positions  .  .  .  . 

Calibration  Positions  .  .  .  . 

Instantaneous  Field  of  View  . 

Total  Cross-Track  Scan  .  .  . 

Total  Scan  Period . 

Integration  Time  (Cross-Track 
and  Calibration  Positions) 

Sync  Mode  . 


Cross-Track  Uadi r 
7 

2-Cold  Space  and  300°K 
14° 

36° 

32  Seconds 

2.7  Seconds 

Auto  Mode  or  On -Sync 


Source: 


Rigone  and  Stogryn  (  1977) 
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Table  2.  Channel  Parameter  Design  Specifications 


Channel 

Polari  zation 

Frequency 

(GHz) 

Bandwidth 

(MHz) 

METD 
( °K) 

1 

50.5 

400 

0.6 

2 

Principal  ly 

53.2 

400 

0.4 

3 

Hori  zontal 

54.35 

400 

0.4 

4 

54.9 

400 

0.4 

5 

Orthogonal 

58.4 

115 

0.5 

6 

to  Channels 

1-4 

58.825 

400 

0.4 

7 

59.4 

250 

0.4 

Source:  Rigone  and  Stogryn  (1977) 

highest  in  the  atmosphere  near  28  km.  The  weighting  functions  result 
from  consideration  of  atmospheric  absorption  due  to  water  vapor  and 
molecular  oxygen  and  indicate  the  approximate  location  in  the  atmos¬ 
phere  from  which  we  assume  most  of  the  energy  in  a  particular  fre¬ 
quency,  that  reaches  the  top  of  the  atmosphere,  originates.  The 
weighting  functions  presented  in  Figure  3  also  have  transmissivity 
corrections  for  antenna  gain  characteristics  included,  as  described 
by  iiipko  (  1979) . 


CHAPTER  3 


ATMOSPHERIC  OPTICAL  PARAMETERS  FOR  SSM/T  FREQUENCIES 

Microwave  frequencies  between  50.0  and  60.0  GHz  are  gener¬ 
ally  considerably  affected  by  several  different  atmospheric  constitu¬ 
ents.  Gases  in  the  atmosphere  absorb  and  emit  microwave  radiation, 
and  liquid  water,  such  as  cloud  droplets  and  rain,  absorb,  emit  and 
scatter  radiation. 

3. 1  Gaseous  Constituents 

The  primary  gaseous  absorbers  in  the  microwave  region  are  H,,0 
vapor  and  02«  Atmospheric  nitrogen  gradually  becomes  a  more  signifi¬ 
cant  absorber  at  frequencies  greater  than  120  GHz,  but  has  little  in¬ 
fluence  in  the  vicinity  of  60  GHz.  Other  gaseous  constituents  occur 
in  very  low  concentrations  and,  in  general,  their  affects  are  very 
smal 1 . 

The  0 2  absorption  spectrum  in  the  microwave  region  is  due  to  the 
molecule's  magnetic  dipole  moment.  Figure  4  presents  the  absorption 
coefficient  of  dry  air  versus  frequency.  The  absorption  coefficient 
(decibels/km)  shows  a  dramatic  peak,  near  60  GHz  and  a  secondary  maxi¬ 
mum  near  120  GHz.  3oth  peaks  are  attributable  to  molecular  oxygen. 

The  characteristics  of  this  absorption  spectrum  were  first  investigated 
by  Van  Vleck  during  World  War  II  and  were  published  in  1947.  Other 
investigators  have  built  on  his  studies  (a  detailed  biography  was  pre¬ 
sented  by  Rosenblum  (1961)),  culminating  in  a  paper  by  Meeks  and 


n 

Lilley  (1963).  For  this  study,  the  absorption  coefficients  for  09 
were  calculated  using  Meeks  and  Lilley  parameters. 

Figure  5  presents  water  vapor  absorption  (dB/km)  versus  frequency 
in  the  microwave  region.  You  can  see  that  the  values  of  absorption 
increase  gradually  through  the  range  of  frequencies  under  investiga¬ 
tion.  It  is  apparent  that  the  absorption  due  to  water  vapor  is  less 
than  that  of  molecular  oxygen  for  all  the  channels  chosen  for  the 
SSM/T  sounder.  Note,  however,  that  even  though  the  oxygen  coefficient 
is  higher,  the  two  values  are  roughly  equivalent  for  the  window  chan¬ 
nel  (50.50  GHz) . 

Studies  of  water  vapor  absorption  (due  to  an  electric  dipole 
moment),  (see  Figure  5),  again,  have  their  roots  in  World  War 
II  projects.  Initial  studies  were  published  by  Van  Vleck  (1947), 
Becker  and  Autler  (1946),  and  King,  et  al.  (1947).  Rosenblum's  arti¬ 
cle  gives  an  excellent  summary  of  work  done  up  to  1961.  HgO  vapor 
absorption  coefficients  were  calculated  using  the  method  of  Barrett 
and  Chung  (1962).  Computer  programs  for  determination  of  both  of 
these  coefficients  were  kindly  provided  by  Falcone  (personal  communi¬ 
cation)  . 

3.2  Cloud  Droplets  and  Rain 

Mie  (1903)  found  an  exact  solution  for  the  scattering  and  absorp¬ 
tion  due  to  homogeneous  spherical  particles,  given  the  complex  index 
of  refraction,  wavelength  of  the  incident  radiation  and  the  particle 
size.  Although  larger  rain  drops  are  known  to  deviate  somewhat  from 
the  spherical  shape,  Mie's  results  still  enjoy  widespread  Dopularity 
and  use,  and  remain  probably  the  best  currently  available  theoretical 
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Figure  5.  Mater  vapor  absorption. 
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tool  for  analyzing  the  effects  of  cloud  droplets  and  rain  in  the 

atmosphere.  The  Mie  results  have  been  adapted  for  computer  use  and 

the  program  used  for  these  calculations  is  a  modification  for  a 

Univac  system  of  the  program  developed  by  Liou  and  Hansen  (1971). 

Outputs  of  the  program  are  s  ,  the  extinction  coefficient  (knf"*)»  w  » 

e  v 

the  single  scattering  albedo,  and  the  coefficients  of  the  Legendre 
polynomial  approximation  to  the  phase  function.  It  is  important  to 
recall  that  the  extinction  coefficient  and  single  scattering  albedo 
calculated  by  Mie's  theory  are  due  to  condensed  water  only.  These 
values  are  then  modified  to  include  effects  of  the  gaseous  atmosphere 
within  the  cloud  prior  to  using  them  in  the  radiative  transfer  equa¬ 
tions.  Tables  of  the  Lengendre  polynomial  phase  function  coefficients, 
extinction  coefficients  and  single  scattering  albedos  for  each  SSM/T 
channel  were  presented  by  Hipko  (1979). 

3.2.1  Indices  of  Refraction 

Savage  (1976)  determined  the  complex  indices  of  refraction  (based 
on  a  least  squares  fit  to  the  data  of  Saxton  and  Lane  (  1962)  given  by 
Hoi  linger  (  1973))  for  pure  liquid  water  at  temperatures  of  +10°  and 
-10°  C  for  several  frequencies  (see  Table  3).  The  indices  of  refrac¬ 
tion  for  SSM/T  channels  (Table  4)  were  derived  from  Table  3  by  linear 
interpolation  (Nipko,  1979).  0°  C  values  were  used  exclusively  in 

this  study. 

3. 3  Drop  Size  Distributions 

Mie  theory  predicts  the  absorption  and  scattering  cross  sections 
for  a  single  spnerical  particle.  Various  investi gators  have  general¬ 
ized  these  results  to  those  for  a  distribution  of  particles  by 
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Table  3-  Refractive  Index (m=mr-im1.)  of  Pure  Liquid  Water 


Frequency 

TS10°C 

T=0°C 

T=-10°C 

(GHz) 

mr  mi 

mr  mi 

mr  mi 

19.35 

5.87 

2.96 

5.20 

2.94 

4.64 

2.81 

37.00 

4.32 

2.60 

3.83 

2.37 

3.46 

2.12 

50.30 

3.76 

2.29 

3.37 

2.03 

3.08 

1.78 

89.50 

3.00 

1.67 

2.76 

1.42 

2.60 

1.20 

100.00 

2.89 

1.56 

2.68 

1.31 

2.54 

1.10 

118.00 

2.75 

1.39 

2.58 

1.15 

2.46 

0.96 

130.00 

2.68 

1.30 

2.53 

1.07 

2.43 

0.89 

183.00 

2.49 

1.00 

2.39 

0.81 

2.33 

0.66 

231.00 

2.40 

0.82 

2.33 

0.66 

2.29 

0.53 

After  Savage  (1976) 


Table  4.  Refractive  Index  of  Pure  Liquid  Water  (SSM/T) 


Channel 

Frequencies 

(GHz) 

50.50 

53.20 

54.35 

54.90 

58.40 

58.82 

59.40 

mr  3. 37 

3.32 

3.31 

3.30 

3.24 

3.24 

3.23 

m.  2.03 

1.98 

1.97 

1.96 

1.90 

1.90 

1.89 

After  Mipko  (1979) 
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integrating  over  the  distribution.  This  area  per  volume  relationship 
has  dimensions  of  inverse  length  and  is  called  the  scattering  (or 
absorption  of  extinction)  coefficient.  Two  models  have  been  selected 
to  represent  clouds  for  the  purpose  of  comparing  with  previous  studies. 
Experimental  data  presented  by  Mason  (1971)  after  results  given  by 
Dunn  (1948)  were  used  to  provide  a  drop  size  distribution  for  non¬ 
precipitating  clouds.  The  data  used  were  for  a  cumulus  congestus 
cloud.  The  second  representation  chosen  was  a  cloud  model  developed 
by  Marshall  and  Palmer  (1948).  This  drop  size  distribution  is  written 

N(r)  =  0.16  exp(-2rS)  (3-1) 

-  21 

where  <5  =  41  R  '  and  R  denotes  the  rainfall  rate  in  mm/hr.  Gunnard 
and  Marshall  (1958),  and  Sekhon  and  Srivastaval  (1970)  have  verified 
Marshall  and  Palmer's  exponential  expression  based  on  raindrops  and 
snowflakes  collected  during  field  experiments. 

3.4  Atmospheric  Profile  and  Cloud  Models  Used 

The  atmospheric  profile  used  in  this  study  was  derived  from  the 
northern  hemisphere  mid-latitude  Spring/Fall  climatology  contained  in 
the  U.S.  Standard  Atmospheric  Supplements,  1966.  The  temperature  and 
height  values  were  interpolated  to  40  standard  pressure  levels. 
Saturation  mixing  ratio  values  for  levels  below  250  mb  were  obtained 
from  a  Skew-T  diagram.  Mixing  ratios  above  250  mb  were  assumed  to  be 
a  constant  3  parts  per  million.  The  Spring/Fall  profile  is  presented 
in  Figure  6. 

A  measure  of  the  liquid  water  content  in  an  atmosphere  was 
obtained  by  the  following  means: 
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Figure  6.  Northern  hemisphere  mid-latitude  Spring/Fall  profile. 
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1)  for  spherical  particles  the  liquid  droplet  volume  per  unit 
volume  of  atmosphere  is 

A  00  *3 

!/*■=•  it  fr  N(r)dr  .  (3-2) 

J  0 

Now  we  may  express  the  liquid  water  content  per  uiit  volume  of  atmos¬ 
phere  as 

W  =  4  N ( r) dr  gm  cm  8  (3-3) 

15  0 

where  r  is  in  cm,  N(r)  is  in  cm'4,  and  the  density  of  liquid  water  is 

_3 

1  gm  cm  .  W  may  be  expressed  for  the  two  cloud  models  as 

W  »  1.16678  x  10'2  gm  cm‘2  km'1  (3-4) 

for  the  Diermendjian  L-Model  and 

W  =  1.77883  x  10'2  •  R'84  gm  cm"2  km'1  (3-5) 

for  the  Marshall-Palmer  model.  Table  5  presents  the  total  mass  of 

2 

liquid  water  within  a  column  of  1  cm  cross  sectional  area  for  the 
cloud  models  used  in  the  microwave  comparison. 
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-2  -2 

Table  5.  Total  Liquid  Water  Content  (10  gm  cm  ). 


Thickness  (km) 

1 

2 

3 

4 

5 

Dei  rmendjian: 

1.16678 

2.33356 

3.50034 

4.66712 

5.83390 

Marshall-Palmer: 

R(nm/hr) 

1 

0.88941 

1.77882 

2.66823 

3.55764 

4.44705 

2 

1.59209 

3.18418 

4.77627 

6.36835 

7.96044 

3 

2.23812 

4.47625 

6.71437 

8.95249 

11.1906 

4 

2.84992 

5.69984 

8.54975 

11.3997 

14.2496 

5 

3.43745 

6.87491 

10.3124 

13.7498 

17.1873 

10 

6.15321 

12.3064 

18.4596 

24.6129 

30.7661 

15 

8.65005 

17.3001 

25.9502 

34.0020 

43.2503 

20 

11.0146 

22.0291 

33.0437 

44.0583 

55.0728 

25 

13.2853 

26.5706 

39.8559 

53.1412 

66.4266 

30 

15.4840 

30.9681 

46.5421 

61.9361 

77.4201 

CHAPTER  4 


BRIEF  REVIEW  OF  MICRCMAVE  RADIATIVE  TRANSFER 
IN  A  PLANE  PARALLEL  , MEDIUM 


4. 1  Clear  Column  Calculations 

The  plane-parallel  model  assumes  that  the  atmosphere  is  made  up 
of  many  layers  of  infinite  horizontal  extent  that  are  homogeneous  and 
in  local  thermodynamic  equilibrium.  This  atmosphere  is  necessarily 
axially  symmetric  about  nadi r-zeni th.  Also,  for  clear  colurm  bright¬ 
ness  temperature  calculations,  there  is  no  source  of  illumination 
at  the  top  of  the  atmosphere. 

The  basic  equation  of  transfer  for  an  absorbing  and  emitting  me- 
di  urn  may  be  wri  tten  as 


dl 

\> 

k  ods 

v 


(4-1) 


where  Iy  is  the  intensity  at  frequency  v,  ky  is  the  absorption  co¬ 
efficient,  p  is  density,  is  the  source  function,  and  S  is  the  slant 
path.  We  change  variables  from  ds  to  dz  such  that  ds=dz/y,  where  z  is 
height,  y  =  cose,  and  9  is  the  zenith  angle.  The  formal  solution  of 
Eq.  (4-1)  is  given  by 

-f  kvpdz'  z  /  kvodz‘ 

I  (z)  =  I  (0)e  0  +  /  J  (z')e  0  kvpdz'  .  (4-2) 

'  0 

where  I  (0)  is  the  radiation  source  from  the  earth's  surface.  For  con- 

v' 

venience,  we  define  the  optical  depth  and  transmittance,  respectively. 
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as  follows: 

t(z' , z)  =  /  k  o  dz,  (4-3) 

z*  J 

T.jtz*  ,z;y)  ■  exp[-T(z'  ,z)/u].  (4-4) 

These  definitions  allow  us  to  rewrite  Eq.  (4-2)  as 

z 

I  ( z)  -  1,(0)  T  (0,z;u)  +  f  J  (z')dT  (z'.z;;.)  (4-5) 

v  v  v  J 

We  may  write  Eq.  (4-5)  to  represent  the  radiance  reaching  the  top  of 
tne  atmosphere  as 

1  -  I  (0)  T  (0,-;u)  +  /  JJz)dT  (z',»)  (4-6) 

J  M  \)  z  -0  ^  ^ 

The  source  function,  J  ,  for  an  atmosphere  in  local  thermodynamic 
equilibrium  may  be  represented  by  Planck's  function.  In  the  frequency 
domain,  it  is  given  by 

3  (T)  =  2hv3/[c2(ehv/kT-l)]  (4-7) 

'J 

where  h  is  Planck’s  constant,  c  is  the  speed  of  light,  and  k  is 
Boltzman's  constant.  In  the  microwave  domain  (hv/KT<<l),  the  Rayleigh 
Jeans  approximation  to  Planck's  function  is  customarily  used  and 
therefore 


exp(hv/kT)=  1  +  hv/kT  . 


(4-8) 


Consequently,  Eq.  (4-7)  reduces  to 
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BV(T)  =  (^-)  T  =  c  T  .  (4-9) 

c 

That  is,  given  a  frequency,  the  Planck  function  is  directly  proportional 
to  the  temperature.  Moreover,  analogous  to  Eq.  (4-9)  the  equivalent 
brightness  temperature  may  also  be  defined  from  the  upwelling  radiance 
in  the  form 

‘v  -  <=vTb(v)  .  (4-10) 

In  a  clear  atmosphere,  radiance  emitted  at  a  level  is  isotropic 
and  makes  two  contributions  to  the  radiant  energy  reaching  the  top  of 
the  atmosphere.  First,  energy  emitted  in  an  upward  direction  travels 
directly  to  the  top  of  the  atmosphere,  suffering  a  Beer's  law  attenua¬ 
tion  due  to  the  absorption  by  atmospheric  gases  along  its  path.  Sec¬ 
ond,  energy  emitted  in  tne  downward  direction  is  attenuated  along  its 
path  to  the  surface.  There,  a  percentage  is  reflected  back  upward 
and  is  again  attenuated  to  the  top  of  the  atmospnere.  Obviously, 
then,  the  contribution  from  the  first  description  to  intensity  reacn- 
ing  the  top  of  the  atmospnere  from  height  z  is 

T^v.z)  =  T(z)Tv(2,»;y)  (4-11) 

The  latter  contribution  may  be  expressed  as 

T  g(v,z)  =  T(z)(l-ev)T;(0,-;u)Tv(0,z;;:)  (4-12) 

d  f* 

Note  Tg  is  the  direct  component  of  brightness  temperature,  Tg,  the 
reflected,  T(z),  the  temperature  at  height  z,  and  T  (z,»;u),  the 
transmission  from  z  to  »  along  path  u.  e  is  the  emissivity  of  the 
surface.  The  atmospheric  contribution  to  microwave  radiance  reaching 


22 


the  top  of  the  atmosphere  may  be  gotten  by  integrating  Eqs.  (4-11)  and 
(4-12)  and  adding  them  together  as 

Tr(v)  =  /  T(z) dT  (z,°°;y)  +  ( 1-e  )T  (0,»;u) 

B  2=0  V  V  J 

00 

•  /  T(z)dT  (0,z;u)  •  (4-13) 

z=0  0 

There  is  one  more  contribution  we  must  take  into  account.  The 
surface  of  the  Earth  emits  in  the  microwave  region  as  a  grey  body,  and 
the  expression  for  this  emission,  attenuated  to  the  top  of  the  atmos¬ 
phere  is  written 

T^(v)  =  £vTsTv(0»^u),  (4-14) 

where  Ts  is  the  surface  temperature  at  z=Q.  The  total  brightness 
temperature  at  the  top  of  the  atmosphere  may  be  written  then  as  the 
sum  of  Eqs.  (4-13)  and  (4-14)  as 

T,(v)  =  f  T(z)dT  (z,»;u)  +  (1-s  )T  ,(0,»;u)  T(z)  dT  (0,z;J 

+  VsV3’^)  •  (4_15) 

Eq.  (4-13)  yields  the  clear  colurtr  brightness  temperature  at  the  top 
of  the  atmosphere.  If  the  transmissivities  have  been  adjusted  for 
antenna  gain  characterise  cs ,  the  Eq.  (4-15)  will  give  us  wnat  are 
commonly  called  "antenna  temperatures" .  Mote  that  we  have  carried 
along  the  dependence  of  the  transmission  on  =  cos  H  (-  =  zenith 
angle)  tnrougnout  this  discussion.  However,  in  clear  column  discus¬ 
sions,  •;  is  normally  assumed  to  be  zero  degrees.  We  have  kept  the 
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nominal  dependence  for  convenience  in  discussing  the  formulation  of 
the  clear  atmosphere  boundary  conditions  for  the  plane-parallel  and 
cylinder  cloud  models. 


4.2  Plane  Parallel  Model  for  Cloudy  Atmospheres 

The  basic  equation  of  transfer  for  a  plane  parallel  cloud  layer 
in  local  thermodynami c  equilibrium  can  be  written  in  the  form: 


di  ( t  j u)  a  i 

— ^ -  =  Iv(t,u)  -  r  /  Pv(u,u‘)  Iv(T,u')du' 

-  (l-iv)  Bv(T(t))  (4-16) 


where  1^  represents  the  monochromati  c  radiance,  p  the  cosine  of  the 
emergent  angle  with  respect  to  tiie  zenith,  x  the  monochromatic  opti¬ 
cal  depth,  P  the  normalized  axially  symmetric  scattering  phase 
function  and  the  single  scattering  albedo.  Bv[T(t)]  is  the  Planck 
function  as  defined  by  Eq.  (4-7).  The  phase  function  can  be  expressed 
as  a  truncated  series  of  Legendre  polynomials.  A  Gauss-Legendre 
quadrature  is  used  to  replace  the  integration  in  Eq.  (4-16).  This 
allows  us  to  derive  a  sat  of  first  order  inhomogeneous  differential 
equations.  seek  homogeneous  and  particular  solutions  of  the  dif¬ 
ferential  equations  as  outlined  by  Chandrasekhar  (1950),  and  the 
complete  solution  of  the  scattered  intensity  from  an  isothermal  cloud 
(temperature,  T  )  for  a  given  discrete  stream,  is  expressed  as 


(4-17) 


where  m  is  the  number  of  discrete  streams  employed,  and  k^  repre¬ 
sent  eigenfunctions  and  eigenvalues  for  the  differential  equations 
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and  the  L  are  a  set  of  constants  determined  from  the  radiation 
m 

boundary  conditions.  Moting  that  the  brightness  temperature  in  the 
microwave  range  is  simply  the  radiant  intensity  divided  by  a  constant, 
Eq.  (4-16)  can  also  be  written  in  terms  of  brightness  temperature. 

The  solution  is  then  given  by 


Tg(v,x;y.j) 


I  a*P(-kmT>  * 

m 


(4-18) 


where  the  constants  L^,  are  different  from  those  in  Eq.  (4-17),  having 
been  obtained  using  brightness  temperature  boundary  conditions.  Tc  is 
the  cloud  layer  temperature  and  is  independent  of  frequency,  v. 

At  the  cloud  top,  the  downward  brightness  is  equal  to  the 
brightness  contributions  from  every  point  in  the  atmosphere  above  the 
cloud  top  attenuated  by  Beer's  Law  down  to  the  cloud  top.  This  can 
be  expressed  as 


;,zt'”’Ji ^ 


/  T(z)dT  (z,zt;-uj 
z-zk 


(4-19) 


where  is  the  height  at  the  cloud  top  and  the  negative  sign  on  u. 
simply  indicates  downward  transfer. 

vlithin  the  cloud  layer,  where  scattering  occurs,  continuity  of 
brightness  values  from  all  directions  is  required.  Thus 


T3 (  .  ,Ui )  -  Tg(v,  r ,  ^ )  l-l  ,2,  .1-1 


(4-20) 


wnere  N  is  the  total  number  of  sublayers  with  the  cloud. 

At  the  lower  boundary  of  the  cloud,  three  brightness  contribu¬ 
tions  are  immediately  apparent: 
a)  the  surface  contribution , 
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b)  the  direct  atmospheric  contribution  from  below  the  cloud,  and 

c)  the  reflected  atmospheric  contribution  from  below  the  cloud. 

A  fourth,  and  perhaps  less  obvious,  brightness  contribution  at  the 
lower  boundary  must  be  considered  for  the  cloud  of  infinite  horizontal 
extent  (Nipko,  1979).  Since  the  reflection  of  downward  brightness  by 
the  Earth's  surface  is  significant  for  microwave  radiation,  the  emer¬ 
gent  brightness  at  the  cloud  bottom  w  11  contribute  to  the  lower 
boundary  condition.  That  is  to  say,  the  solution  to  the  radiative 
transfer  through  the  cloud  affects  the  boundary  conditions  used  to 
obtain  the  solution.  Nipko  (  1979)  suggested  an  iterative  approach  to 
correct  the  lower  boundary  condition.  If  we  include  this  additional 


term,  we  may  express  the  lower  boundary  condition  as 


VwvV  *  e»TtV°’W  ♦  zi0  T(2)dTv(2’2b;-i 


)  V°'zb 


4b 

)[  '  T(z)dT.  (0,z;-u, ) ! 
z-0  ''  1 


T3 (  J ’ 2 1 ; - u  1 '  Vyi '  1  u!'  3,Zb’ui^ 


(4-21) 


where  T  (u..)  is  defined  to  be  the  top/down  throughput  of  the  cloud  for 
the  stream  defined  by  q.  and  z^  is  the  height  of  the  cloud  base. 


i0(-i)  =  T3  ( z  b  ’  ”  u  i  >  /T3  ( v  ;  u  i ) 


(4-22) 


The  solution  of  the  microwave  transfer  equation  given  by 
Eq.  (4-13)  is  applicable  only  to  isothermal  and  homogeneous  cloud 
layers.  To  aoply  the  transfer  solution  to  thick  clouds  it  is  neces¬ 
sary  to  divide  tne  cloud  into  sublayers  for  which  the  isothermal  and 


26 

homogeneous  assumptions  are  reasonable.  By  matching  the  brightness 
temperature  continuity  equation  for  brightness  temperature  in  each 
sublayer,  a  set  of  linear  equations  with  unknown  coefficients  can  be 
obtained  which  can  be  solved  by  standard  matrix  inversion  techniques. 
Similar  procedures  have  been  used  by  Liou  (1975)  to  evaluate  the 
transfer  of  solar  radiation  in  inhomogeneous  atmospheres  and  by 
Feddes  and  Liou  (  1977)  to  investigate  the  transfer  of  spectral 
infrared  radiation  in  cloudy  atmospheres. 


CHAPTER  5 


RADIATIVE  TRANSFER  IN  A  CYLINDRICAL  MEDIUM  - 
A  SPHERICAL  HARMONICS  APPROACH 


Radiative  transfer  in  cylindrical  geometry  will  be  investigated 
for  a  one  dimensional  cylindrically  symmetric  infinite  cylinder  and  a 
two  dimensional  infinite  cylinder.  The  source  function  was  chosen  so 
as  to  represent  terrestrial  radiation.  Note  that  the  frequency  depend¬ 
ence  index,  v,  is  omitted  for  simplicity. 


5.1  Basic  Equations  and  Simplifying  Assumptions 

Monochromatic  radiative  transfer  can  be  described  by  the  equation 


dl  _ 


SQds 

e 


I  -  J 


(5-1) 


I  is  the  monochromatic  intensity  (or  radiance),  J  is  the  source  func 
tion,  3g  is  the  extinction  coefficient,  and  s  is  a  distance  in  space 
(Chandrasekhar,  1950).  Since 


dl  _  * 
ds  ' 


(5-2) 


Eq .  (5.1)  becomes 


1_  9I(s  ,5;t) 
c  3t 


I(s,D:t)  -  J(s 


t) 


(5-3) 


where  T  is  a  unit  vector  specifying  the  direction  of  scattering 
through  a  position  vector,  s.  We  further  reduce  the  equation  by 


28 


assuming  a  steady  state,  homogeneous  medium  and  we  may  now  write 

-  0>7)I(U)  =  1(1,3)  -  J(U).  (5-4) 

4e 

In  the  infrared  and  microwave  spectral  regions,  the  source  function 
J  for  a  scattering,  emitting  and  absorbing  atmosphere  may  be  written 
as 

J(s,n)  =  £  /  I(s,o)P(o,5-)dn'  +  (l-S  )3v(T)  (5-5) 

4m 

where  5,  is  the  single  scattering  albedo,  B^(T)  is  the  Planck  function 
of  temperature,  T,  at  wavenumber  v,  and  P  is  the  scattering  phase 
function  which  is  normalized  such  that 


4tt 


/ 

a-rr 


p(n,^* )dn‘  =  i 


(5-6) 


The  source  function  terms  represent  multiple  scattering  and  emission, 
respectively ,  from  the  assumed  homogeneous  medium,  and  the  parameters 
P,  and  3g  are  assumed  to  be  independent  of  the  position  vector,  s. 
We  now  write  Eq.  (5-4)  in  cylindrical  coordinates  as 


3_ 

3z 


]i(s,n) 


I(S,5)  -  J(s.3)  . 


(5.7) 

The  relevant  aeometry  of  the  cylinder  is  illustrated  in  Fig.  7. 
r  is  a  position  vector  from  the  center  of  the  cylinder  to  a  point  on 
the  cylinder  wall,  2  is  described  by  9,  the  angle  from  the  z  axis, 
and  p,  the  difference  between  the  projection  of  ,i  onto  the  plane 
containing  r,  and  ,  the  angle  of  the  polar  coordinates  of  r.  As  can 
be  seen  from  Fin.  7,  the  components  of  2  may  be  written  as 
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3  =  sin0  cos-P 
sine  si  no 

*r  "  " 

3  =  cosS  (5-8) 


We  now  make  the  further  assumption  that  the  z  axis  of  the  cylinder  is 
infinite  and  the  problem  becomes  independent  of  z.  Eq.  (5-7)  becomes 
(substituting  in  Eqs.  (5-8)) 


37  Jr  -  57 


-]  I (r ,3)  =  I (r ,3) 


/  Kr.^Ptf.sWdn'  +  O-ujB  (T) 

4tt  0  J 


(5-9) 


We  will  seek  approximate  solutions  for  Eq.  (5-9)  using  the  SDherical 
harmonics  method.  We  expand  the  phase  function  as 


P(3,3') 


»  l 

1  1  Zj 

7=0  m=-7 


Ym* 

l 


(3‘)  Y*  (3) 


(5-10) 


The  assumed  form  for  the  solution  for  the  intensity  may  also  be 
written  as 

00  £ 

I(r,3)  =  l  l  Im(F)Y^(3)  +  B  (T).  (5-11) 

7=0  m--l  c  J 


where  *  denotes  the  complex  conjugate, 
sion,  and  the  spherical  harmonic  Ym(3) 


Tj?  are  coefficients  of  expan- 
is  defined  as 

(cose)eim*  .  (5-12) 


~he  Pm(cos0)  are  the  associated  Legendre  polynomials  described  by 
Janke-Emde  (1945).  These  spherical  harmonics  satisfy  the  orthonor¬ 
malization  relation  such  that 
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/  y™(h)y^*  (n)dn  = 


~Z  -m 
2a+l  Ca  °3 


[5-13) 


0  m 

where  5  and  o'  are  the  Kronecker  delta  functions. 

a  3 

useful  relationships  are  the  addition  theorum 


Some  additional 


l 


P^cose)  =  l  YrJ(Q)  Y'J  (fl) 
m=-Jl 


(5-14) 


and  the  recursion  formula 


sinee1'^  Y™(Q)  =  (  /U-m-1 )  U-m)  Y1^  (3) 


m+1 


m+1 


-  S(l+m+\)(l+m+2)  Y'''+‘(n)} 


(5-15) 


and 


sin9e"1<5  Y^(5) 


2J+]  {  /U-m+1TU-m+2y  Y^;j  (3) 


-  /U+m-D(J.+m)  Y™"](n)}  .  (5-16) 

First,  we  will  investigate  Eq.(5-9)  for  the  case  of  cylindrical 
symmetry,  i.e.,  the  intensity  is  a  function  of  r,  6  and  and  then 
will  return  to  the  more  general  problem  where  the  intensity  is  a 
function  of  r(r,p^),  9  and  <p. 

5.2  One-Dimensional  Model  (Cylindrical  Symmetry) 

5.2.1  First  Order  Solution 

We  rewrite  Eq.  (5-9)  to  reflect  cylindrical  symmetry  as 


(Kofink,  1959) 
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d3  _  (a+g-l)1/2  (a+g)1/2 
a  2(2a-l) 


( 5-20d ) 


Some  details  of  this  procedure  may  be  found  in  the  Appendix.  By  noting 
that  t  =  3gr  we  may  rewrite  Eq.  (5-19),  for  simplicity,  as 


a3C^- 

ordx 


- 


"Sir  -  w  - 


csri 

ordt 


^-]I3+](t) 

t  a-1 


+  o3[4-  - 

aLdx 


— HS“!(t)  =  v  I3(i)  -  (1-5  )B  (T) 
t  -  a-1  a  a  '  v  v  v 


(5-21) 


where 


V  =  1 

a 


03  ,0) 

v  a 
2a+l 


(5-22) 


For  a  first  approximation  to  the  set  of  partial  differential 
equations,  (5-15),  we  select  the  zero  and  first  modes  and  set  the 
1-2,  ...<*>  modes  to  zero.  Allowing  these  degrees  of  freedom  we  obtain 
the  four  basic  equations: 


A°0[£*i]l’<T)  - 

-  n-s;)Bv(T) 
(5-23) 

-i  d  To,  v  _  T-i,  ,  Vi  T-i,  , 

C1  dt  X1  ^  3  *1  ^ 

(5-24) 

r°(x)  -  o 

(5-25) 

_1  d  ,0/  >  Tl,  s  Vi  Tl/  i 

D1  dT  *0^  "  ‘  3 

(5-26) 

Solving  for  i]  and  in  Eqs.  (5-24)  and  (5-26),  respectively,  and 
substituting  into  Eq.  (5-23),  and  taking  the  homogeneous  portion  only, 
yields 
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2 

r  ^  2 
dx 


+  T 


dx 


+  x"n2Ig(x)  =  0 


(5-27) 


which  is  Bessel's  equation  of  order  zero  with  a  parameter  n.  Here 


(5-28a) 


where 


❖i +  “oV 


a 


•JO  JO, 

V  1 


( 5 -28b ) 


The  solution  to  Bessel's  equation  may  be  written  in  terms  of  the  modi¬ 
fied  Bessel  functions  as 


Tg ( x )  =  ClI0(^)  +  C2k0(I)  (5-29) 

where  •  =  IQ  and  kQ  are  the  zeroth  order  modified  Bessel  functions 
of  the  first  and  second  kind,  respectively,  and  c-j  and  Cg  are  constants 
to  be  determined. 

There  are  two  boundary  conditions  to  be  applied  for  the  first 
order  aDproximation  with  cylindrical  symmetry.  First  we  recognize 
that  the  intensity  is  finite  at  -  =  0.  Since  kg(— )  approaches  minus 
infinity  at  x  =  0,  we  must  have  c^  =  0.  Therefore 

Io(x)  =  C]I0(X)  •  (5'30) 

The  second  condition  is  on  the  outer  boundary  and  will  be  discussed  in 
general  in  the  next  section. 


5.2.2  Extension  to  Higher  Order  Approximations 

Vie  assume,  as  oroposed  by  Kofink  (1959),  that  the  system  of 
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partial  differential  equations  denoted  in  Eq.  (5-21)  is  satisfied  by  a 
linear  combination  of  modified  Bessel  functions 


where  c^  ,  are  arbitrary  constants  and  P™  and  Q™  are  eigenfunctions 
of  the  eigenvalue,  X .  Once  again,  applying  the  boundary  condition  at 
t  =  0,  we  find  ~  0  and 

=clxO^Vf)  (5-32) 

Through  the  recursion  properties  of  the  modified  Bessel  functions 
we  find  that  it  is  possible  to  transform  Eq.  (5.21)  into  recursion 
relations  for  the  eigenfunctions  P™(,\).  The  recursion  formula  are: 

I_m(x)  =  Im(x)  (5.33a) 

‘  ^+1^^  =  7  lm^'1  (5.33b) 

1  dIn/^ 

^Wx>  +  Wx>}  -HST-  •  (S-33C) 

From  Eqs.  (5.33)  we  find 


Therefore  we  may  write 
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.iripin+l  ( .  •,  Rmpm-1  r .  < 

ArW'J  B2P2+11^ 


-  CV+I(  )  +  ) 

L  A,*  X  C  ».  -  1 


.  ?*'  1 


{5-35} 

By  eliminating  the  square  roots  in  the  coefficients  by  the  formula 


P*(‘)  =  (22+1 


(-D 


U-m)/2 


.  ^U-m) !  U+iTT  pm,  n  ,, 

:  * ;  v  /  /  >  v  3 


5-36 ) 


the  recursion  formula  for  the  P^1  ( • )  is  transformed  into  a  simpler 


relationship  in  R^(') 


(2+m+l  )R^I  ('■ )  +  +  (■'--m)R^!(  )  *  U+m)Rm_1  ' 


r  ; 


=  2(22+1)  y,  R.  (  , ) 


'5-37 ) 


The  solution  for  Eq.  (5.37)  may  be  written  (according  to  Kofink) 

(5-38) 


r”(.)  =  j”(z)  g<2)(A) 


where  the  suoerscriot  (z)  represents  the  main  solution  (z=Q)  or  the 
zth  complementary  solution  (z  =  1,°=).  We  will  discuss  the  main  solu¬ 
tion  first. 

The  main  solution  satisfies  the  recursion  formula 
U+Dg^jt-)  +  ig^}(;-)  =  (22+i)v,;.g(?)(  )  ,  (5-39) 

and  3^°^  Hi.  The  set  of  gi°'  (  )  begins  with  l  =  3  with 
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(0), 


arbitrarily  set  to  1.  We  write  out  the  first  few  g :  '()  in  terms  of 
qg°V-)  as 
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40>t) 

9<°»(  ) 


g<°>(> 


=  1 

=  *r0  gg°V- ) 

/3  .2  1  x  (0) ,  > 

(2  '  '  2^g0  ^ 

r5  .3  ,5  .  2  \  i  (0),.  , 

C2  VlV  '  '6  3  0^  g0 


(5-40) 


It  is  convenient  to  write  Eqs.  (5-40)  in  matrix  form  as 
1  0  0  0... 


I  -V 

I 

;  1 


0 

0 


2  -5y2-  3 


0 

0 

4 


.(0) 


(0) 


i  t 


^1 

g<01 

i 


=  0. 


L  • 


(5-41) 


"nis  system  of  equations  for  the  main  solution  can  be  carried  out  to 
the  order  of  aoDroxiration  desired.  A  nontrivial  solution  to  Eq.  (5-41) 
r.ay  be  obtained  by  setting  the  determinant  equal  to  zero.  Vie  then 
find  the  roots  of  the  characterise  polynomial,  wnich  are  the  set  of 
.  There  will  be  (n  +  1 ) / 2  pairs  of  roots  for  the  main  solution  since 
n  is  normally  taken  to  be  odd.  The  roots  of  the  characteri Stic  poly¬ 
nomial  may  be  gotten  directly  for  n  =  1,3.  We  find 


for  n  =  1 


and  the  ‘‘our  roots 
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where  a  =  105  vgf ’  b  =  "(27  YqY  +  35  y^y.  +  28  y^y^)  and 
c  =  9. 

The  characteristic  polynomial  for  higner  order  n  was  determined 
by  the  following  scheme.  The  tridiagonal  matrix  is  of  the  form 


(5-42) 

and  the  characteristic  function  is  given  by  the  sequence 


1 

al 

a2fl  *  b2b3f0 
a3f2  ‘  b4b5fl 


j 

t 


b2n-lb2n-2fn-2 


For  n-5,  we  are  no  longer  able  to  solve  directly  for  the  eignevalues, 
and  therefore  have  employed  Newton’s  iteration  scheme  to  numerically 
calculate  the  roots. 
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Mark  (1945)  showed  that  one  additional  degree  of  freedom  in  the 
solution  is  created  every  time  n  is  increased  by  two.  These  "comple¬ 
mentary"  solutions  are  reflected  in  Eq.  (5-38)  and  the  recursion  for¬ 
mula  for  the  complementary  solutions  corresponding  to  that  for  the  main 
solution  given  by  Eq.  5-39  (Kofink,  1959)  may  be  written 

U+l-2z)g^}[42)]  +  (Z+2z)g =  (2i+l );. <z)y;g<2)  , 


?,  =  2n 


n  =  1 ,  2 ,  3 


[5-43) 


The  coefficients,  £^,  of  Eq.  (5-38)  may  be  determined,  then,  by  the 
recursion  relations 

U+m+l)2^]  ^  =  2(Z-2n+1)SlJ  ^  (5-44a) 


(:.-■)£]  (z)  ♦  ( l+m )  (z)  =  2(i+2n)3™  (z) 


[5.44b) 


For  example,  beginning  at  ( T ,m )  =  (0,0)  with  gQ  arbitrarily  set  to  1, 

0(1 ) 

we  may  add  the  arbitrary  constant  R,,'  to  the  main  solution  at  (i,m)  = 


q  ['£0)]  =  g^0)L  £0)>R°(1)!:  jn)]  • 


(5.45) 


’  1 1  l  i;  1 )  1  *  the  first  complementary  solution,  begins  a  new  series  of 
functions  with  l  =  2,  that  obey  the  recursion  formula  in 
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Eqs.  (5-43)  and  (5-44).  Starting  with  arbitrarily  set  to  1, 

we  solve  the  recursion  formula  (Eq.  (4-43))  straightforwardly  until  we 
reach  the  order  of  the  approximation  we  desire  (£  =  n).  For  an  example, 
we  write  out  the  matrices  needed  for  a  fifth  order  solution  (n=5). 

For  the  main  solution  we  have 


-Y 


.(0) 

G  '  k 

1 

0 

0 

Q 

0 


1 

_3y  x(0) 
JY1  k 


2 

0 

0 

0 


0 

2 


-5yoa 


•  (0) 
2Ak 


3 

0 

0 


0 

0 

3 


7v  '  W 

"7f3Ak 


4 

0 


0 

0 

0 

4 


-9y 


(0) 


4  k 

5 


-11 Y 


0 

0 

0 

0 

5 

„  ,(0) 


5  k 


,(0)  ! 
0  ! 
,(0)  ; 


g2 


(0) 


J  L 


.(0) 

’3 

,(0) 


40)j 

(5-46) 


=0 


For  the  first  complementary  solution  we  have 

■  ( 1 ) 

'2 

Q(1)l! 

"3  ,1 

nO)  ! 


-5v 

Y  j 

5 

1 

0 

0 

r 

j 

3  J 

6 

2 

0 

0 

-9v -1  ^ 

3 

i 

0 

0 

4  J 

7 

-nv,)n 

95 


(1) 


For  the  second  complementary  solution  we  have 


-9y 


(2) 


4'  i 
9 


-llY5"i 


(2) 


q(2 ) 

! 

94 

n  ( 2 ) 

g5 

=  0. 


=  0  .  (5-47) 


(5-48) 
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Once  we  have  determined  the  ,\£z\  we  may  calculate  the  functions 

and  the  modified  Bessel  functions  I _[t/.\£z^]  (please  see 
^  k>  m  k 

Eq.  (5-32)).  We  now  turn  our  attention  to  the  calculation  of  the 
unknown  coefficients,  c^[\^]  by  means  of  applying  the  boundary  con¬ 
dition  at  the  outer  surface. 

5.2.3  Boundary  Conditions  at  the  Outer  Surface 

From  the  preceding  discussion,  we  have  seen  that  there  are 
(n  +  1  )/2  pairs  of  eigenvalues  in  the  main  solution  and  (n  -  l)/8 
pairs  of  additional  roots  occuring  in  the  complementary  solutions.  We 
add  these  to  find  (n  +  l)(n  +  3)/8  pairs  of  involved.  Therefore, 
we  look  for  the  same  number  of  outer  boundary  conditions  to  determine 
the  unknown  constants 

Historically,  the  rigorous  boundary  condition  first  applied  to 
this  type  of  problem  was  (in  the  neutron  transport  literature)  the  so 
called  “Vacuum11  boundary  condition,  e.g.  ,  no  energy  incident  on  the 
boundary  or 

I (~  ,u)  =  0.  (5-49) 

It  is  obvious  that  this  condition  cannot  be  met  exactly  by  any  finite 
expansion  of  the  intensity.  Mark  (1945)  suggested  that 

I(t,ui )  =  0  (5-50) 

where  the  intensity  was  expanded  in  a  truncated  series  of  Legendre 
polynomials  and  the  u.  were  the  >./ 2  positive  roots  of  P  (c.)  =  0. 
Marshak  (1947)  suggested  setting 

1 

j  I  (  "  ,  M  )  Po  '  4_1 
0 


(;-)dc  *  o 


(5-51) 
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where  the  intensity  is  expanded  in  the  sane  truncated  Legendre  form  to 
approximate  the  ’Vacuum"  boundary.  It  is  worthy  of  note  that  Marshak's 
proposal  seems  to  be  much  more  widely  applied  than  Mark's  for  a  zero 
incident  intensity  distribution  (e.g.,  Pomraning,  1964;  Gelbard,  1968; 
Dave  and  Canosa,  1974). 

In  this  study  we  are  interested  in  the  situation  when  there  is  an 
intensity  distribution  incident  on  the  outer  boundary 

Ur.S)  =  f  (£2) .  (5-52) 

In  this  situation,  a  generalization  of  Marshak's  integral  conditions 
appear  to  enjoy  popular  acceptance.  These  conditions  may  be  written 
1 

/  P2n+1(u)  (Kt.u)  -  f (u) }  dp  =  0.  (5-53) 

Since  the  incoming  intensity  distribution  incident  on  the  boundary, 
f(p),  is  presumed  known,  one  may  carry  out  the  various  angular  inte¬ 
grations  . 

For  this  study,  an  analogous  method  was  used  to  approximate  the 
boundary  conditions.  Rewriting  Eq.  (5-10),  for  convenience,  and  re¬ 
calling  that  we  are  dealing,  at  present,  with  the  cylindrically 
symmetric  case,  and  recognizing  the  form  of  I™(t)  (Eq.  (5-32)),  we  have 


(5-54) 


I(t,3)  »  B  (T)  +  j  I 
Z-0  m =-Z 


where 
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vT(n)  +  C-l )m  Y"m(5) ,  m/0 

yJ(o)  m=0 


(5-55) 


Note  that  we  have  combined  the  spherical  harmonics  using  the  relation¬ 
ship  =  (-1)m  (Kofink,  1959).  Now  we  choose  a  weight 

function  (as  in  Marshak's  method),  multiply  both  sides  of  Eq.  (5-54)  by 
it,  and  integrate  over  the  appropriate  domain,  The  choice  of  the  weight 
function  should  be  related  to  the  spherical  harmonics  and  is  taken  to 

be  Y^(9,o)  where  (a+3)  =  even.  The  domain  of  the  incident  hemisphere 

a 

is  0  *  6  *  tt  and  j  ,  Carrying  out  the  indicated  operations, 

we  may  write 

m  3  it/2  n  „ 

/  /  B  (T) Y^(9 ,p) 7^(9 ,c|))sin9dedd) 

0  VI 


where  the  right  side  of  the  equation  is  the  known  incident  intensity 
distribution.  The  evaluation  of  the  integrals  of  the  produces  of  the 
spherical  harmonics  becomes  quite  tedious  for  higher  order  approxima¬ 
tions.  For  this  study,  the  integrals  were  accomplished  analytically 
up  to  and  including  n=5.  A  general  numerical  scheme  for  these  inte¬ 
grations  has  been  developed  by  Ou  (1979,  private  communication). 
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Continuing  with  n=5,  for  example,  we  look  for  (n+l)(n+3)/8  or  6  outer 
boundary  conditions  for  our  one- dimensional  study.  These  will  be 
obtained  by  accomplishing  the  operation  described  by  Eq.  (5-56)  with 
the  following  weight  functions 

yj(fi)  =  y](&)  -  Y^(sJ)  =  -  Sz  pj(u)  cosq 
=  Y^(ft)  -  Y"1^)  =  -  J|p^(u)  coso 
ip)  -  -  P33(S)  cosaq 

YjCiT)  =  -  j^r  Pg(u)  cosq 

Y5^  =  "  /50W  P5^  cos3<1, 

Y5^  =  "  j 90 7, 2 00  P5^  cos5t?>  (5-57) 

5.3  Two  Dimensional  Model 
5. 3. 1  General  Solution 

We  recall  Eqs.  (5-9) -(5-15) .  In  the  one-dimensional  case  we 
simplified  I(r,£)  to  I(r,*).  For  this  instance,  we  take  I(r,^)  to  be 
I(r,4>r,3),  where  r  is  the  polar  coordinates  of  any  "slice"  across  the 
assumed  infinite  cylinder.  'We  follow  essentially  the  same  procedure 
as  before  by  substituting  in  the  expansions  for  the  intensity  and  the 
phase  function,  and  utilizing  the  orthonormalization  property.  We 
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then  carry  out  the  operation 

/  Eq.  (5-9)  •  Y3*(H) dfi  .  (5-58) 

4rr  * 


Once  again,  using  the  recursion  formula,  Eqs.  (5-15)  and  (5-16),  we 
obtain  the  new  set  of  partial  differential  equations. 

3  r_d_  ,  (3+1)  1  d  ■,  T  3+ 1  t  ,  \ 

(t,4>r) 

(  t  ,  d>r) 

(T,i>r) 

(5-59) 

Details  of  this  integration  may  be  found  in  Appendix  A. 


j.  ui  i  I  l  uy  i  r 

P 


R3  r±  .  (3-1)  ,  1  d  ,  t3-1 
a  at  t  it  dij>  ct+1 


r3  ;_d_  ,  ( 3+1 )  _  A  ^  i 

at  t  "it  dq>r "  a-1 


,  d3  i-A  -  ie^Li  ♦  J.  <L,  ,3-1 

a  ldx  t  it  d<}>  a-1 

r 


=  Y  I3  (t,4>  )  -  (1-5  )B  (T) 
'a  a  v  r  ^  j  \t  7 


5.3.2  First  Order  Approximation 

Selecting  the  zero  and  first  order  modes  we  may  write  from 
Eq.  (5-59)  (homogeneous  portion  only): 


A 


0 

0 


J __d_i  t  1  / _  .  j  R0 
it  dq  '  1  i(T>l?r)  -  Bq 
r 


(t,  *>r)  = 


(5-60) 
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1  r  d  _ 1_  d  ,  T  0  /  ,  \ 

'C1  lddr  "  ix  d<$r  '  !0 


I^1(T,<Dr)  [1 


y1]  (5-61) 


*!-o 


(5-62) 


°]  { df; +  it  d|;}  lilted 


(5-63) 


We  solve  Eqs.  (5-61)  and  (5-63)  for  I^and  I j ,  respectively,  and  sub¬ 
stitute  into  Eq.  (5-60).  After  some  algebraic  manipulation,  and 
noting  that 

0  1  .  0  -1  „  1 
A0  D1  "  B0  C1  '  6  * 


we  may  write  our  new  expression  as 

1  r  d  ,  1  d  ,  1^  d  -t  .  0  /  ,  \ 

177  7^]I°(T,V 


Vo(T'V 


(5-64) 


We  now  see  that  Eq.  (5-64)  is  equivalent  to 


K,72lQ(x,br)  -  K2Iq  =  0  , 


This  is  a  time-independent  diffusion  equation  expressed  in  cylindri¬ 
cal  coordinates,  where 


“  3y, 


K-2  =  Tg 


Separation  of  variables  is  used  to  solve  Eq.  (5-64).  We  let 


Ig(x,dr)  =  R(r)t>  (i»r) 


( 5-65) 
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and  substituting  Eq.  (5-65)  into  Eq.  (5-64)  and  separating  variables, 
we  find 


t2  d2R  ,  t  dR  2  K2  2 
I  ^+R  dt  '  T 


( 5-66) 


and 


]_  d2i> 
T 


Eq.  (5-67)  may  be  written  as 


=  n 


h  ♦  n2  t  .  0 


and  the  solution  is 


(5-67) 


5  =  cos(nor)  +  Bn  s i n ( n d> r) 


(5-68) 


Eq.  (5-66)  may  be  expressed  as 


.2  d2R  ,  ,  dR  ,2  K2  2n  „  _  0 

:  — *•  +  T  -3—  -  [7  rr~  +  n  )  R  -  0 

o7  dT  Ki 


K, 


This  is  the  modified  Bessel  equation  with  a  parameter  l  ='7—. 

1  '  1 
before,  we  let  l  =  — ,  and  the  solution  to  Eq.  (5-66)  is 


As 


R  =  C1>n  In(y)  +  C2>r)  Kn(j),  n  unrestricted.  (5.69) 
From  Eq.  (5-65)  the  complete  solution  is 


I0<T-=r)  *  J0  *  C2  1[An  cos(Mr>  *  Bn  s1n<"V5 

(5-70) 

Two  boundary  conditions  are  applied  at  this  time  for  simplicity.  As 
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before,  we  must  have  C-,  =0  since  ;<  ->■  as  r-*0.  Next,  since  ‘(0)  = 
■(2-t)  on  any  "slice"  of  the  cylinder,  we  must  have  n  be  an  integer. 
We  may  now  rewri  te  Eq.  (5-70)  (combining  and  renaming  the  constants) 
as 

-  0° 

IqU,^)  =  cos^n<j!P  +  Bn  sin(%)]  •  (5-71) 

Eqs.  (5-61)  and  (5-63)  may  now  be  solved  for  I^(t,or)  and  lj(r,or), 
respectively.  We  recall  the  first  order  approximation  expansion 
for  the  intensity  (homogeneous  portion  of  Eq.  (5-10))  may  be 
expressed  as 


I(t,.'.)  =  Ig(T,ir)Yg(0,o>)  +  l](T,©r)Yj(3,b)  +  l/d^Y^e.b) 

(5-72) 


We  may  write  I ^ 1  ( r , and  l](r,^r)  as 


/[  y  i  n 


l  {2X  [In-l(T)  +  In+l(l)^An  cos(ndr) 


+  Bn  sin  ( nd>r)  ]  +  I  n(j)  j  sin(n©r)  +  Bn  cos(n$r)]} 

(5-73) 


a-  d 


■  jr~  l  <2T  +  Wi^ 

vc  Y|  n 

+  8n  sin(n©r)]  -  In(^)  j  [-A^  sin(npr) 


+  3  cos^ni  )  ]  . 

n  t 


(5-74) 


We  substitute  Eqs.  (5-71),  (5-73)  and  (5-74)  into  Eq.  (5-72  and  may 
write 
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l(r,<Pr,d,<b)  =  l  In(y)  CAn  cos(nqr)  +  sin(nd>r 


&  Yl  n 


I  =«(n»r) 


y|> 

it 


n[-An  sin(n<?r)  +  Bn  cos(n<fr)]}  Vj(Q) 


&  y-j  n 


l  <lt  »„.,(T>  +  I„+1Cf>^A„  cos(nV) 


vi 

it 


n  [-/^  sin(n<ir)  +  Bn  cos(nd>r)]}  Y^(3; 


To  combine  terms  we  add  and  subtract  y](3)  and  Y^(? 

yj(3)  -  y"1^)  =  -  Si  p](u)  coso 


and 


Y J ( 3)  +  Y"  1  (5)  =  -  S2  i  P -1  ( u )  si 


,-l 


1V“'  '  '1 


sing 


where  P^(u)  are  the  associated  Legendre  functions. 
Eq.  (5-75)  to  reflect  these  combinations, we  have 

I(t,i>r,9,o)  =  l  cos(nor)  +  8n  sin(nqf) 


*  ■—  nh  t'n-ltf)  +  cos("ir)  * 

vcY|  n 

•  CYj(T)  -  Y-Vn)]; 


)] 

+  Bn  sin(no^)  ] 

+  3n  sin(n$r)] 
(5-75) 

’)  as 

(5-76a) 

( 5  -  75b) 
Rewriting 

J 

3n  sin(ntr)] 
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i  Ur) 

+  -  I  {  — *r~n[-An  sin(nd)r)  +  Bn  cos(rv&r)] 

v?  Y i  n 

•  CYj(^)  +  Y^(n)]>  •  (5-77) 

Upon  substituting  Eqs.  (5-76)  into  Eq.  (5-77),  we  arrive  at 

I(t,4  ,0,$)  =  T  In(y)  [An  cos(n<})r)  +  Bn  sin(nj>r)] 
n 

P I ( u) COSi 

-  -W -  “s<nV  +  sn  s1n(nV; 

1  n 

p](u)sin<6 

-  WfO&n  c°s<n9r>  -  \  s1n<n»r)]- 

I  n 

(5-78) 

■Note  that  we  have  made  use  of  the  recursion  properties  of  the  modi¬ 
fied  Bessel  functions  (Eq.  (5-33b))  in  the  last  term  of  Eq.  (5-78). 

At  this  point,  the  unknown  constants,  An  and  Bn>  may  be  determined 
by  applying  the  outer  boundary  conditions.  The  discussion  of  these 
boundary  conditions  will  be  held  in  abeyance  until  we  have  presented 
the  extension  of  the  two-dimensional  model  to  higher  orders. 

5.3.3  Extension  to  Higher  Order  Approximations 

To  generate  higher  order  approximations  for  the  two  dimensional 

problem,  we  look  for  a  general  solution  for  the  ).  From  the 

0  1  -1 

form  of  the  first  order  parameters,  Ig(r»Pr)>  I-|(~**r)  and  I-]  (t,or), 
and  by  comparison  with  the  one  dimensional  problem,  we  are  led  to 


assume  the  form 


Substituting  this  solution  into  the  homogeneous  form  of  Eq .  (5-59) 
and,  again,  noting  the  recursion  properties  of  the  modified  Bessel 
functions,  we  may  transform  Eq.  (5-59)  into  a  recursion  relation  for 
the  eigenfunctions  as  before  (see  Appendix).  '/Je  find 

that  we  can  form  the  same  recursion  formula  for  the  as  in 

the  one  dimensional  case  -  that  is,  the  are  invariant  for 

the  one  dimensional  and  two  dimensional  Droblems.  We  may,  therefore, 
determine  the  eigenvalues  '-[z  *  and  the  eigenfunctions  Pm[Vz']  by  the 
same  methods  presented  previously. 

5.3.4  Boundary  Conditions  at  the  Outer  Surface 

As  was  written  previously,  the  form  of  the  approximation  to  the 
intensity  is 

to  ' 

Kt,')  =  '  7  Im(-)Ym(C)  +  B.  (T) 

;=0  m=-r.  ''  " 

Substituting  in  tne  expression  *ar  lT(t),  Eq.  (5-79),  we  have 
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®  i 


kt.3)  •  i  y  i  p"[->t(2)]  y  <[ - *- 


’n+jTTzP  +  ’n-m^TTT1 


1=0  m=-l  . (z) 


*  [An[X^]  cos(n<?r)  +  Bn[X^]  sin(na)r)] 


In+m  '  rn-m[^TFT] 

+  1  [ - - - 2 - ][Bn[X^2)]  cos(nor) 


-  An[X^]  sin(n<tr)]}  Yj(5)  +  BV(T) 


(5-80) 


For  example,  the  expansion  for  the  third  order  approximation  is 
straightforward  but  tedious  and  is  written 


I(T,n)  =  B.  (T) 


■'(z)  n 

'  k 


l  InC7(zTjCAn[^2)]  cos(n®r)  +  Bn[:^z)]  sin(ncr) 


•  p§L^z)iv°a)  ♦  p]cf  )]c 


In+l[rriT]  +  In-l['TT7] 


[An['^z)]  cos(npr)  +  Bn[A^z)]  sin(n®r )][y|  (5)  -  Y^tf)] 


+  pj[^z)][ 


W:T7T]  " 

'k 


Xn-1  C7fzT] 

- - - ][Bn['^z)]  cos(n<Pr) 


-  An[:^z)]  sin(nor)][i*|tf)  -  Y'1  (X;,.]  +  P°[^z)]  InC4jy] 


[y42^  cos(%)  "  Bn['Jz]]  sin  (n;  }] 


r '-L  ? 


ran)] 


'"*2^  *  ^-zCrrry^ 
fy-rll - * - , - A—ira  r,U) 


][A  J]  cos  (  n:  ) 


nL  k  -1 


sin(ncr)][Y22(n)  +  y2^)]  +  p2^)-. 


In+2^.CzT^  ~  In-2t~iiy- 


-][Sn[^z)]  cos (no  )  -  A  [ ^ 2 '  ] s i n ( n ;■ 


r  i  r v2 


-  r/ag  *  p  ['  z>][ - 5 - , _ t _ ■ 


(,n[:k''j  «»(%)  *  ?„[]2)]  Sln(n.-  )i[y'(C)  .  y;1,'); 


r  3 ' 


W:^  *  ^-iLrrrr1 

3L  k  jC  5  ][3n[‘]-'j  cos(n;  ) 


V'kz,3  »i»(«r):c«y'e)  *  ?;’(";)] .  oft  <*>: 


T  I" 

‘n+B1- 


3L.,'U)'  +  rn-3^TTT' 


~  ■^'V'k2^  cos'n';r) 


+  Bn['.  2)]sin(n: 


Cy?(3) 


-  V(fi)]  +  p^[xfz)][ - A 


n+3[7fiT^  -  ^-3^] 


LV\  ^  cos(n<J>r)  -  AnE^Z^]  sin(npr)]|j(y2(o)  +  y;3(Jg;; 


W  • 


^ote  that  we  have  used  the  relationship 

f^2)]  =  <-om 


(5-81) 


to  simplify  the  expansion. 

-q.  (5-31)  may  be  further  simplified  by  collecting 


common  factors 


BV(T)  + 


*  L(z)  |  C\C' k*’]  cos("Pr)  +  Bn[i<2)]  sin(%)] 

'k 

n  ,  .  In*lki(2)J  +  'n-jf-pyJ 

'  '■  ‘nCyiT1  P0['k  1  V0®  f  p  [>■  Z’][ - ^ 


•  tr|(3)  .  yj'&j  *  $3)  ♦  p|r.(*>3 


■n^rtTyl  T  'n-2^7z}l 

[ - 5 - 5 - —  ICfftk*)  -  C2(:)]  *  p1[.<z)] 
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(5-82) 


As  in  the  one  dimensional  model,  the  outer  boundary  condition 


on  this  problem  is  that  the  "moments  of  flux"  across  the  boundary 
are  continuous.  We,  therefore,  are  led  to  multiply  Eq.  (5-S2)  by  the 
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same  weight  functions  as  before  and  integrate  over  the  domain  of  the 
incident  hemisphere.  The  right  hand  side  of  Eq.  (5-82)  will  be  the 
assumed  known  incident  intensity  distribution.  This  intensity  distri¬ 
bution  may  be  calculated  by  assuming  that  the  infinite  cylinder  is 
embedded  in  a  non-scattering  (clear  atmosphere)  plane  parallel  medium. 
See  Chapter  4  for  details  of  clear  column  intensity  calculations. 

A  double  Gauss  quadrature,  utilizing  16  streams  of  intensity  in  each 
angle  (9  and  $  mapped  onto  the  interval  (-1,  1)),  was  used  to  inte¬ 
grate  the  right  hand  side  of  the  equation  after  it  was  multipled  by 
the  appropriate  weight  function.  Note  that  multiplying  the  right  hand 
side  by  the  weight  function  Y^(6,i>)  and  integrating  over  the  domain 
results  in  what  we  normally  call  the  flux  -  but  we  will  term  it  the 
"first  moment"  of  flux.  Y^*?  ,<?)  is  the  "second  moment"  of  flux,  and 
so  on,  to  avoid  confusion.  Calculations  of  the  moments  of  flux  were 
performed  at  18  equally  spaced  intervals  on  the  cylinder  wall  between 
0  and  ~r.  Recognizing  that  the  incident  moments  of  flux  are  periodic 
in  9  ,  and  noting  the  form  of  the  left  hand  side  of  Eq.  (5-82),  we 
are  lead  to  express  the  right  hand  side  moments  in  a  Fourier  series. 

In  fact,  the  right  hand  side  moments  are  even  functions  on  (0,  ~)  and 
may  be  reduced  to  a  set  of  cosine  series.  The  right  hand  side  moments, 
therefore,  may  be  written  as 

Right 

Hand  =  (dQ  +  £  dp  cos(nctr) 

Slde  n=1  y](9,j);  ^(0,0);  etc.  .  . 

The  coefficients,  d^  and  dn,  of  the  cosine  series  of  each  set  of  13 
moments  were  determined  numerically  (Powers,  1972)  by  the  algorithms 
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d0  =  i4f (x0)+f (xi )+'  ’ ,+f (Vi  )4f (xs) 5 


and 


dn  =  (xQ)+f(x1  )cos— x1  +•••+  Jr  f(x$)(cos 


m 

a 


xs) , 


n=l ,  •  ■  •  ,  s  -1 , 

where  the  interval  (0,rr)  is  divided  into  s  equal  subintervals  with 
endpoints  0=Xq,  X-j  ,  •••,  Xs  =  it.  Returning  our  attention  to  the 
left  hand  side  of  Eq.  (5-82),  we  note  that  when  we  multiply  by  the 
weight  function  and  integrate  over  the  incoming  hemisphere  (<$>  from 
tt/2,  3tt/2  and  0  from  0,tt),  about  half  of  the  terms  involve  integra¬ 
tions  of  the  form 

3tt/2 

/  sind  costpdq  £  0. 
tt/2 

In  fact,  these  are  the  terms  with  the  coefficient 
{  BnC:'{Zh  cos(nor)  -  An[X^z)]  sin(ncr)}  . 

For  an  example,  we  may  now  express  the  first  moment  of  flux  as 

3tt/2  1T  i  zr?T 

/  /  Eq.  (5-32).  V  (9,d)  sin  dedd  =  S  (T)  Y  Y' 

tt/2  0  v  i  u 


+  l  )  [A  I'S ^  cos(nd  )  +  B  [ ]  sin(n<D  )] 
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where 

/  2  TT 

=  */  I  [Y "(j?)  +  (-1  )mY:m(Q)2[Y®(^)  +  (-1  )mY"^)]sined9d43 
tt/2  0 

Thus  we  see  that  Eq.  (5-82)  is  a  set  of  X^  simultaneous  equations 
for  any  given  n.  We  may  determine  the  unknown  constants  by  use  of 
the  orthogonality  of  the  Fourier  series.  That  is 

Ao[xk2)]  =  27  /  (5’83)  ddY 

-7T 

and 

Bn[A][z)]  =  l  f  Eq.  (5-83)  sin(n^f)  dop 

and 

An[4z)]  *  \  f  Eq.  (5-83)  cos(n*r)di>r  . 

-7T 


(5-83) 


CHAPTER  6 


COMPUTATIONAL  RESULTS 

6.1  Comparison  with  Cubic  and  Plane-Parallel 

Cloud  Models  in  the  Infrared 

Calculations  were  performed  using  the  two-dimensional  infinite 
cylinder  model.  Mie  theory  was  used  to  generate  the  single-scattering 
parameters  for  a  typical  water  cloud  at  the  chosen  wavelength  of  10  urn. 
Values  employed  were  a  single-scattering  albedo,  3  ,  of  .63  and  coeffi¬ 
cients,  ,  of  the  expansion  of  the  phase  function.  The  are  approxi¬ 
mately  2.3,  3.4,  3.6,  3.4  and  2.7  for  £=1,5.  Note  that  absorption  due 
to  water  vapor  was  not  included  in  the  infrared  calculations. 

Comparisons  were  made  with  two  different  models.  Liou  (1973)  gave 
results  for  a  plane-parallel  cloud  using  the  16  stream  discrete  ordinate 
method.  Liou  and  Ou  (1979)  presented  calculations  for  radiative  trans¬ 
fer  in  Cartesian  coordinates  for  cubic  and  rectangular  clouds.  A 
spherical  harmonics  expansion  of  the  scattering  phase  function  and  in¬ 
tensity  was  used  for  the  finite  clouds  and  the  expansion  was  truncated 
at  the  first  mode  (e.g.,  £=0,1)  and  may  be  called  a  first-order  solu¬ 
tion. 

Physically,  both  cloud  models  were  assumed  to  be  isothermal  with 
temperature  T.  and  were  situated  above  a  horizontally  infinite  diffusely 
emitting  plane  whose  temperature  is  T  .  Note  that  no  atmosphere  was 
present.  The  cylinder  cloud,  for  comparison  purposes,  was  configured 
in  the  same  manner.  Boundary  conditions  for  the  cylinder  cloud  were 
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calculated  at  a  series  of  equally  spaced  points  on  a  cylinder  cross 

section  as  described  in  Chapter  5.  Because  the  assumed  black  surface 

is  the  only  source  of  radiation  outside  of  the  cloud,  we  see  that  tne 

downward  flux  density  will  be  zero  at  the  cloud  top  point  (or=-r),  the 

upward  flux  density  will  be  ttB.  (T  )  at  the  cloud  bottom  point  (<?  =0), 

and  the  flux  density  arriving  on  the  circular  side  of  the  cloud  will 

be  a  smooth  function  with  values  between  0  and  mB  (T  ).  Note  that  the 

j  s 

cylinder  cloud  is  situated  in  the  atmosphere  so  that  the  polar  angle, 

6  ,  equals  zero  at  the  nadir  direction  and  it  at  the  zenith  direction. 

Emergent  intensity  from  the  cloud  models  is  a  quantity  of  some 
theoretical  interest.  Liou  and  Ou  (1979)  presented  average  intensity 
values  emerging  from  the  top  and  bottom  faces  of  cubic  and  rectangular 
clouds.  These  values  were  compared  to  emergent  values  from  the  plane- 
parallel  cloud  (Liou,  1973).  The  cloud  top  point  and  cloud  bottom 
point  are  the  only  points  on  the  cylinder  cross  section  that  lie  in 
the  same  plane  as  the  top  and  bottom  "faces"  of  the  other  two  models. 
Hence,  the  emergent  values  at  those  points  were  chosen  for  comaprison 
purposes.  Figure  3  details  the  intensity  (per  units  of  -B  (T  ))  from 
the  cloud  top  and  bottom  as  a  function  of  emergent  angle  for  vertical 
cloud  optical  depths  of  1  and  10.  Note  that  the  fifth  order  intensity 
values  are  presented  here.  In  each  case  the  finite  cloud  models  pro¬ 
duced  straight  up  values  and  straight  down  values  that  were  larger  and 
smaller,  respectively,  than  the  plane-parallel  model.  For  upward  emer¬ 
gent  angles  near  90°,  each  model  predicted  a  smaller  value  than  the 
plane-paral lei  cloud. 

Remote  sensing  applications  demand  an  average  intensity  value 
weighted  over  an  area  that  is  "seen"  by  an  observer  or  satellite  from 
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Figure  8.  Incident  intensity  versus  emergent  angle  per  -B  (T  )  in 
the  upward  and  downward  direction  at  10  :j.  : 
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a  far  distance.  Physically,  then,  a  weighted  average  intensity  would 
be  the  observed  value.  Because  of  the  scanning  method  of  infrared  and 
microwave  instruments,  clouds  may  be  viewed  at  several  angles  of  obser¬ 
vation.  For  example,  a  cloud  viewed  at  45°  obviously  will  present  a 
different  surface  than  a  cloud  viewed  at  nadir.  For  this  purpose  we 
may  define  a  weighted  average 

N 

Iw(r,4>r)  =  l  Ii(x,4»r;(fri)  •  w.($.) 

where  Iw  is  the  weighted  average  intensity,  I.  are  the  N  emergent  in¬ 
tensity  values  in  direction  $  ,  the  angle  of  observation,  and  oi  are 
the  equally  spaced  angles  at  which  the  emergent  intensity  values  are 
chosen.  Liou  and  Ou  (1979)  also  define  an  area-weighted  upwelling  in¬ 
tensity  emergent  from  the  surface  of  a  cubic  cloud.  Comparison  of  re¬ 
sults  of  weignted  average  intensity  calculations  for  different  angles 
of  observation  are  presented  in  Figure  9.  Again,  note  that  fifth  order 
intensity  solutions  for  the  cylinder  cloud  are  graphed  for  comparisons. 
The  figures  for  optical  depths  1  and  10  show  that  the  weighted  average 
intensity  viewed  at  nadir  for  the  cylinder  cloud  is  much  closer  to  the 
plane-parallel  result  than  to  the  cubic  cloud.  The  large  value  por¬ 
trayed  at  nadir  for  the  cubic  cloud  may  be  partially  explained  by  re¬ 
calling  that  the  cubic  cloud  solution  is  only  a  first  o^der  approxima¬ 
tion.  A  first  order  truncation  of  the  cylinder  weighted  average  inten¬ 
sity  also  produces  a  large  value  at  nadir  -  equal  to  .375.  Mote  also 
that  a  third  order  turncation  produces  a  nadir  value  of  .265  for  the 
cylinder  cloud,  smaller  than  the  fifth  order  solution.  From  this  we 
may  observe  the  expected  oscillatory  behavior  of  the  spherical  harmonics 


j  average  intensity  versus  angle  of  observation  at  10  pm 
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solution.  When  the  angle  of  observation  approaches  45°  we  see  that  all 
of  the  models  predict  nearly  the  same  value.  Wear  90°  both  the  cylin¬ 
der  and  cubic  clouds  exhibit  similar  behavior.  As  the  optical  deptn  of 
the  cloud  is  increased,  the  multiple  scattering  and  emission  effects 
within  the  cloud  manifest  themselves  in  increased  energy  emerging  from 
the  sides  of  the  cloud.  Obviously  these  results  must  be  tempered  by 
the  fact  that  we  have  not  considered  an  intervening  atmosphere  between 
the  cloud  and  the  sensor.  It  does  appear,  however,  that  the  relative 
observation  angle  and  the  cloud  position  within  a  sensor  scan  may  be 
important  in  the  remote  sensing  of  cloud  composition  from  satellite. 

Another  parameter  of  interest  in  infrared  radiative  transfer  is 
cloud  emissivity.  Cloud  emissivity  may  be  investigated  by  removing  tne 
surface  emission  from  the  transfer  equation.  Effectively  this  means 
that  the  intensity  incident  on  the  boundaries  of  the  cloud  models  is 
zero.  Figure  1 0  shows  the  average  emissivity,  defined  as  t.ie  ratio  of 
the  average  upward  flux  density  to  -rB.  (T  ),  as  a  function  of  vertical 
optical  thickness.  Emissivity  values  were  calculated  for  the  first, 
third,  and  fifth  order  approximations  to  the  cylinder  cloud  model.  As 
we  saw  in  the  weighted  average  intensity  calculations,  the  emissivity 
curves  for  the  cylinder  cloud  manifest  the  expected  oscillatory  behav¬ 
ior  due  to  the  order  of  the  spherical  harmonic  approximation.  Since  we 
would  expect  the  fifth  order  emissivity  to  be  slightly  higher  than  the 
true  value,  and  because  the  value  of  the  fifth  order  emissivity  ap¬ 
proaches  a  number  slightly  less  than  1,  the  cylinder  model  results  tend 
to  support  the  conjecture  by  Liou  and  Ou  (1979)  that  clouds  of  finite 
extent  snould  not  be  considered  as  black  bodies  even  when  their  optical 
deptr.s  are  large.  Therefore,  the  cloud  top  temperatures  inferred  from 
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iaure  10.  Averaqe  emissivity  versus  vertical  ootical  thickness 
at  10  jn . 
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optically  thick  clouds,  utilizing  the  10  urn  window,  will  oe  too  nign  if 
the  cloud  emissivity  is  not  accounted  for. 

6.2  Cylindrical  Cloud  Model  Configuration 

in  Realistic  Atmospheres 

The  scattering,  absorbing,  and  emitting  cylinder  is  placed  Hori¬ 
zontally  in  an  absorbing  and  emitting  plane-parallel  atmosphere.  The 
surface  of  the  Earth  is  assumed  to  be  a  horizontal  plane  of  infinte  ex¬ 
tent  that  emits  diffusely.  Since  we  have  assumed  an  infinite  cylinder, 
each  cross  section  of  the  cylinder  is  like  any  other  and  we  may  look  at 
one  such  cross  section.  Calculations  for  the  boundary  conditions  were 

oerformed  at  18  points  for  0C:  where  i  =0  is  the  nadir  direction 

r  r 

with  respect  to  the  cloud  and  o  =~  is  the  zenith  direction.  In  general 
one  could  place  a  tangent  plane  on  the  cylinder,  find  streams  of  inten¬ 
sity  at  various  Gauss  latitudes  in  5  and  ;  and  integrate  using  Gauss 
Quadrature  to  find  the  inward  moments  of  flux.  The  problem  is  simpli¬ 
fied  somewhat  by  the  infinite  cylinder  and  Plane  parallel  atmosphere 
assumptions.  As  a  result,  by  taking  into  account  the  position  angle, 

:  ,  we  may  project  onto  the  cross  section  plane  and  find  the  stream 
angle  through  the  plane-parallel  atmosphere  that  we  need.  We  may  now 
proceed  to  calculate  the  intensities  in  those  stream  angles. 

The  intensities  in  each  stream  are  now  calculated  in  much  the  same 
manner  as  in  the  boundary  conditions  for  the  plane  parallel  model  with 
some  minor  difference  due  to  the  position  of  eacn  tangent  to  the  cylin¬ 
der  cross  section.  Figure  11  illustrates  the  geometry  and  it  is  seen 
that  only  the  very  top  point  will  have  no  streams  tnat  contain  surface 
effects,  and  that  only  the  very  bottom  point  will  have  no  streams  due 
to  the  atmosohere  only.  The  incoming  streams  of  intensity  for  the  very 
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w’vio''  include  terms  for  surface  emission,  atmosoheric  emission,  and 
surface  reflection. 
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an  irtans'ty  stream  wit*  tires  components  to  a  stream  wit-  :r’y  one 
component.  ~c  adequately  represent  tnis  gradient,  ‘c  Gauss  '•toes 
c-csen  in  cotn  5  ana  o  angles  resulting  in  ’23  different  ct^eatmc  •'  . 
for  eacn  point.  These  intensity  / alues  were  t.nen  ■'elated  to  eacn  ,  : 
pair  and  after  multiplying  py  tn.e  appropriate  weight  functor  ;as  de¬ 
scribed  before,,  we  integrate  utilizing  a  double  Gauss  cuadrature  to 
get  tne  moments  of  flux. 

. 3  Comparison  wltn  Plane-Parallel  Model 
for  tne  CMS?  Microwave  Frequencies 

Sensitivity  studies  were  performed  for  microwave  frecuencies. 

Cloud  thickness  of  1,  2,  2,  4  and  2  km  were  used  in  simulating  rainfall 
rates  of  2,  3,  4,  IQ,  15,  20,  25  and  20  mm /hr  oy  use  of  the  Marsrall- 
Da’-ner  drop  size  distribution  as  described  in  Chapter  2.  The  s:oud 
oase  height  was  fixed  at  a  height  of  1  km. 

"he  weignted  average  intensity  for  tne  cylinder  model,  as  de- 
scribed  in  tne  infrared  study,  was  used  to  compare  with  plane-oaral lei 
results  with  one  modi fication.  Circe  eacn  eoually  spaced  point  in  a 
real ’ Stic  atmcsonere  on  tne  cylinder  cross  section  is  at  a  different 
optical  depth,  tne  streams  were  attenuated  to  the  too  of  the  atmosphere 
prior  to  applying  the  weighted  average. 

rigures  12,  13  and  14  display  the  variation  of  brightness  tem¬ 
peratures  as  nadir  at  tne  to p  of  tne  atmosphere  over  ’and  ;-.=0.97:  and 
ocean  51’  surfaces  for  frequencies  of  50.50  ,  53.20  and  o4.35  GHz 

’Channels  1-3),  roppecti vel y .  ~he  cloud  temperature  is  assumed  to  be 
222. If'7  ’n  all  cases.  Note  tnat  brightness  temperature  at  the  too  of 
tne  atmcsonere  r"or  a  frequency  of  54.90  GHz  is  affected  almost  negngi- 
,i  /  py  me  thickest  cloud  considered  in  this  study.  41  so  tne  weighting 
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functions  of  frequencies  53.-0,  53.32,  and  59.40  GHz  peak  suf'i-c’ent'y 
high  in  the  atmosphere  trat  their  origntness  temperatures  are  not 
affected  by  cloud  models,  'nerefore,  we  nave  chosen  not  to  include 
these  frequencies  in  the  sensitivity  studies.  For  convenience,  tne 
abbreviations  PP  and  CL  will  be  used  when  referencing  tne  pi  are-para i lei 
and  cylinder  models,  respectively. 

At  5G.5Q  GHz  over  land  (Figure  12),  the  CL  model  displays  an  ini¬ 
tial  rise  in  origntness  temperature  (Tg)  for  very  small  liquid  water 
contents,  followed  by  a  sharp  drop,  and  then  levels  off  to  a  slightly 

decreasing  trend.  The  Tg  finally  reaches  a  value  of  263.8  at  a  liquid 

„2 

water  content  of  .7742  gm  cm  .  This  liquid  water  content  corresponds 
to  tne  most  optically  thick  cloud  modeled  -  a  5  km  cloud  with  a  rain¬ 
fall  rate  of  30  mm/hr.  The  RP  model  responds  similarly  to  the  CL  model 
over  a  land  surface.  There  is  no  initial  rise  in  the  Tg  values,  how¬ 
ever,  and  the  final  Tg  value  reached  is  near  259.35,  a  difference  of 
approximately  4. 

Over  ocean,  the  50.50  GHz  frequency  values  (Figure  12)  for  ooth 
models  experience  a  strong  rise  from  the  clear  column.  The  CL  model, 
as  in  the  land  case,  has  a  higher  value  throughout  tne  spectrum  of 

liquid  water  contents.  The  difference  in  the  two  models  appears  to  be 

_2 

largest  at  .7742  gm  cm  -  a  Tg  difference  of  approximately  3. 

At  53.20  GHz  (Figure  13),  the  Tg  differences  between  the  two  models 
cor  large  liquid  water  contents  is  practically  negligible.  There  are 
some  differences  in  response,  however,  for  liquid  water  contents  of 
less  tnan  .30  gm  cm  .  Over  land,  both  models  show  an  initial  rise 
in  T,  from  the  clear  column  /alue  but  the  CL  model's  rise  is  much 

•j 

larger  'aoout  Tg=8).  Over  water,  the  PP  model  drops  sharply  and  then 
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'-ises  to  become  larger  in  value  than  tne  CL  model  (aoout  T.^1)  near  a 
liquid  water  content  of  .15  gm  cm'6  before  reaching  nearly  tne  same 
value  at  .30  gm  cm  .  The  T3  value  for  tne  CL  model  over  water  rises 
initially  to  a  peak  value  of  near  259.  and  then  drops  slowly  to  tne 
values  of  tne  other  curves. 

At  54.35  GHz  (Figure  14),  the  response  of  the  Tg  to  varying  cloud 
tnicknesses  and  rainfall  rates  is  much  less  than  for  the  50.50  and 
53.20  GHz  frequencies.  Note  that  the  scale  for  the  Tg  has  been  consi¬ 
derably  reduced.  The  behavior  of  this  frequency  is  similar  to  53 .20  GHz 
in  that  the  responses  to  the  various  forcing  functions  reaches  a  near 

equilibrium  value  at  a  liquid  water  content  of  .30  gm  cm'6.  Below 

_2 

.30  gm  cm  ,  however,  we  do  see  some  differences.  The  CL  model  ever 
both  surfaces  experiences  an  initial  rise  from  the  clear  column  value 
of  ~D  of  approximately  7.  The  PP  model  over  land  also  shows  this  ini¬ 
tial  rise  in  Tg  but  on  the  order  of  4.  The  Tg  of  the  PP  model  over 
ocean  gradually  rises  to  tne  equilibrium  value  near  .30  gm  cm  “  liquid 
water  content. 

Note  that  the  onl y  change  in  the  transfer  model  between  the  ocean 
and  land  cases  studies  is  a  change  in  the  surface  emissivity  -  assumed 
for  the  comparison  to  be  .97  for  a  land  surface  and  .51  for  an  ocean 
surface.  This  change  of  surface  emissivity  essentially  affects  two  of 
the  boundary  conditions:  1^  the  surface  emission  and  2)  the  atmospheric 
contribution  reflected  from  the  Earth's  surface,  which  is  approximately 
16  times  greater  over  ocean  than  land.  Note  also  that  tne  50.50  GHz 
frequency  over  ocean  results  for  both  models  are  generally  higher  than 
tne  clear  column  value  and  under  score  the  usefulness  of  this  fre- 
cuency  in  'flagging''  cloudy  conditions  that  may  cause  a  deterioration 


in  the  accuracy  of  temperature  retrievals  cased  on  microwave  scunner 
data . 

5.4.  Interpretation  of  Results  prom 

Cylindrical  Cloud  Model 

Inis  section  will  address  several  different  topics.  First  tne 
cloud  top  brightness  will  be  split  into  its  various  components  in  an 
attempt  to  discern  why  the  CL  model  responds  as  it  does.  Next  tne 
average  intensity  at  the  top  of  the  atmosphere  will  be  treated  from 
the  remote  sensing  point  of  view.  The  third  topic  focuses  on  tne 
importance  of  multiple  scattering  in  the  DMSP  microwave  frequencies. 
Finally  the  spherical  harmonic  solutions  of  different  orders  are 
analyzed  with  zero  boundary  conditions. 

Figure  15  snows  in  a  general  way  what  mechanisms  are  at  work  to 
cause  energy  gain  and  loss  within  a  cylindrical  cloud.  Energy  gains 
are  realized  from  multiple  scattering,  droplet  and  gaseous  emission. 
Energy  is  lost  from  the  discrete  streams  emerging  from  the  cloud  from 
droplet  and  gaseous  absorption  and  single  scattering.  The  interplay 
of  these  mechanisms  determines  the  solutions  emerging  from  tne  cylin¬ 
der.  The  forcing  functions  on  the  problem  are  the  temperature  of  the 
cloud,  through  the  Planck  function,  and  the  boundary  radiation  streams 
that  illuminate  the  cylinder. 

To  detail  how  the  cloud  interacts  witn  the  clear  column  atmosphere, 
we  examine  the  three  contributors  to  cloud  brigntness  temperature:  1} 
emission  by  the  Eartn's  surface,  2)  atmospheric  emission,  and  3)  emis¬ 
sion  from  witnin  the  cloud.  We  may  examine  tne  effect  of  eac.n  of  these 
terms  independently  by  "turning  off"  two  of  tne  three  contributions 
simultaneously.  Trms  may  oe  accomplisned  by  I)  eliminating  tne  surface 


emission  ~or.tr'  cut'  Gn  r>om  tne  ce'ow-tne-nor-i  zon  streams  of  intensity 
'n  toe  ccurdary  conditions,  2)  setting  t.ne  atmosone^ic  cent ri outicr  to 
ze^o,  and  3)  setting  tne  temperature  of  tne  cloud  to  3f<. 

Figures  ISanc  17  present  tne  component  ccntri Duti ons  to  cloud  too 
brightness  for  50. 5G  GHz  for  tne  CL  and  PP  mode’s,  respectively.  Mote 
tnat  all  values  are  fifth  order  solutions.  Results  are  presented  fo r 
a  5  <m  thick  cloud  with  a  temperature  fixed  at  233 . 7  6  '  K .  t_and  (-c=.a7  j 
and  water  (e=-.51)  surface  calculations  are  grapned.  Note  tnat  tne 
trends  for  all  of  the  various  components  to  cloud  top  origntness  are 
roughly  the  same  for  tne  CL  and  PP  models.  Minor  differences  in  tne 
trends,  however,  may  be  used  to  explain  the  responses  of  tne  CL  and  ?° 
models  seen  previously  in  Figure  12.  At  liquid  water  contents  greater 
tnan  .40  gm  cm  ",  the  cloud  top  origntness  and  tne  Tg  decrease  very 
slowly.  The  solution  due  to  the  cloud  alone  is  nearly  constant  in 
both  models.  The  difference  is  seen  in  the  contributions  from  the  sur¬ 
face  and  atmospheric  effects.  The  surface  effect  in  the  3P  model  is 
effectively  zero  at  cloud  liquid  water  contents  greater  than  .4  gm  cm*". 
The  CL  model,  however,  continues  to  snow  an  effect  due  fo  tne  sjrface 
throughout  tne  spectrum  of  cloud  liauid  water  contents.  Simi1arly,  the 
atmospheric  effect  for  the  CL  model  at  large  liquid  water  contents  is 
slightly  larger  tnan  the  PP  model.  These  contributions,  plus  tne  fact 
tnat  tne  CL  cloud  solution  is  approximate! y  3,5K  higher  tnan  tne  PP 
model  explain  the  difference  seen  in  this  region  in  Figure  12.  3elow 
.3  gm  cm  “  cloud  liquid  water  content,  the  atmospheric  and  surface 
effects  continue  to  have  larger  values  for  tne  CL  model.  These  larger 
values  are  oartialiy  offset  since  the  PP  model  cloud  solution  rises 
’icre  quicKly  tnan  the  CL  solution.  The  net  effect,  though,  is  for  tie 
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water  contents,  tne  atmcscnere  arc  s..r~ace  tcrtr-put-tr >  ve 
-rca'  *•/•»'•  a  laro  do  ret  dreo  fast  ercugr  tc  c sneer  sate  *;<•  v- 

se  •  r  one  cloud  solution.  Tnerefpre,  we  soserve  */e 
r*  se  •'■  to  2  value  larger  tr.an  clear  column  *n  t;;$  r®c*cr . 

gu'es  13  and  '3  present  ere  component  ccr.tr:  out* on.s  to  c'cuo 
too  octrees  "or  53 .23  3Hs  *cr  tre  3w  and  -c  models,  '•esoecti  /el  / . 

-3  seer  orev'ously  in  Fig.  5,  tre  solutions  *or  cloud  ’ ’  ;uid  Nats'-  con¬ 
tent:  greater  t.nan  .35  cm  cm*'  are  near'/  -dentical  -‘or  octr  mcce  1  s 

,  _  ,  _n 

over  octr  surfaces,  .-or  /a ; ues  isss  tnan  . 3C  gm  cm  *  t*ere  are  eras'--: 

d-*rerences  in  the  responses  of  tre  models.  In  rig.  1 3  tre  3L  model 

s tews  ar  ;nitiai  r-lse  rrcm  tre  clear  eclumn  values  over  ootr  surfaces, 

and  t-en  a  crco  toward  tre  eouil i  Sri  urn  ocsit'on.  ~rese  rises  are 

seem-rgi/  caused  oy  tne  observation  tnat  tre  slooe  of  t«e  CL  atmos- 

;rer-c  and  sur-acs  solutions  tends  to  ce  s'Hgntly  more  ncrizonta'  as 

tre  o'oud  1-cuid  water  content  is  reduced,  wr.ve  tre  cloud  solution  ; : 

r-smg  most  a  c  i  a  1  y  .  As  cloud  liquid  water  -s  increased  tne  piece  o“ 

fe  atmosoneric  and  surface  graces  increase  and  cause  t^e  net  solution 

to  decree se  s'! owl /. 

As  discussed  and  defined  previously  in  tne  infrared  sensitivity 
study,  an  important  oarameter  'rom  a  remote  sensing  point  of  view  is 
tne  weighted  average  intensity,  tne  intensity  tnat  an  observer  or 
satellite  would  see  f'-om  a  -'ar  distance.  r:gure  3*.  presents  t"e 
we-grtad  -average  ongrtness  temperature  at  tre  too  o*  tre  atmosonere 
as  a  *urct: on  of  t-e  angle  of  observation  “'or  5  <m.  si  suds  of  optical 
deotr  3.5  and  12.  o/er  a  1  and  surface  “iff:  order  solutions  .  Also 
presented  :r  - ;q.  20  are  solutions  for  a  oiaie-oaral lei  cloud  witn  tre 
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same  configuration.  The  cloud  temperature  was  assumed  to  be  265.41°K. 
The  CL  solutions  for  both  cloud  optical  thicknesses  gradually  decline 
as  the  angle  of  observation  increases  to  near  60°  from  nadir.  From 
60°  to  90°,  the  CL  values  drop  more  rapidly,  and  the  solution  from  the 
optically  thinner  cloud,  although  higher  in  value  from  0°  to  60°,  drops 
below  the  solution  of  the  optically  thicker  cloud.  The  PP  cloud  solu¬ 
tions  decline  slowly  to  near  30°  and  then  more  quickly  from  30°  to 
90°.  At  nadir-zenith  we  might  expect,  as  discussed  earlier  in  the  IR 
study,  that  the  fifth  order  solution  of  average  intensity  for  the  CL 
model  might  be  slightly  higher  than  the  true  value  because  of  the 
oscillation  inherent  in  a  spherical  harmonics  solution.  The  PP  solu¬ 
tion  is  16th  order  and  should  be  reasonably  accurate  for  an  infinite 
slab.  Between  60°  and  90°  the  solution  for  the  optically  thick  CL 
cloud  declines  less  rapidly  than  that  of  the  optically  thin  cloud.  The 
net  result  is  that  at  90°  the  thick  CL  cloud  has  a  higher  solution. 

The  explanation  may  be  found  in  Fig.  21.  The  CL  values  presented  here 
are  intensities  at  the  view  angle,  averaged  as  in  the  IR  case.  That 
is,  the  intensity  values  are  averaged  without  attenuation  to  the  top 
of  the  atmosphere  and  therefore  represent  an  average  intensity  leaving 
the  cloud.  It  is  apparent  that  the  optically  thick  CL  cloud  has  an 
average  intensity  leaving  the  cloud  that  increases  as  the  angle  of 
observation  approaches  90°.  The  optically  thin  CL  cloud  solution  de¬ 
clines  as  the  view  angle  approaches  90°.  Apparently,  as  the  CL  cloud 
optical  depth  is  increased,  the  multiple  scattering  and  emission 
effects  within  the  cloud  begin  to  manifest  themselves  in  increasing 
amounts  of  energy  emerging  from  the  sides  of  the  CL  cloud.  This  effect 
cannot  be  seen  in  the  PP  model  because  of  its  geometrical  limitation. 
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Even  though  the  fifth  order  solution  from  the  optically  thick  CL 
cloud  is  rising  near  90° ,  the  corresponding  brightness  temperature  at 
the  top  of  the  atmosphere  is  decreasing.  The  reason  may  be  that  the 
CL  cloud  is  immersed  in  a  clear  plane-parallel  atmosphere.  Therefore, 
as  the  angle  of  observation  approaches  90°,  the  amount  of  atmosphere 
slant  path  between  the  observer  and  the  cloud  approaches  infinity. 
Therefore,  the  amount  of  energy  reaching  the  observer  from  the  cloud 
should  approach  zero.  This  effect  is  shown  in  Fig.  22  for  a  5  km  CL 
cloud  of  optical  depth  12.  Obviously  for  frequencies  that  fall  at  or 
near  the  absorption  band  center,  the  atmospheric  optical  depth  between 
the  cloud  and  the  sensor  is  already  so  large  that  little  or  no  energy 
from  the  cloud  reaches  the  sensor.  In  "wing"  and  "window"  frequencies, 
however,  relative  observation  angle  and  cloud  position  within  the  sen¬ 
sor  scan  may  be  of  importance  in  the  remote  sensing  of  cloud  composi¬ 
tion  from  satellite. 


Discussion  of  the  importance  of  multiple  scattering  in  the  DMSP 
microwave  frequency  region  is  of  paramount  importance  because  if  there 
is  little  or  no  multiple  scattering  occurring,  the  solution  of  the  PP 
model  is  greatly  simplified.  The  equation  of  transfer  in  a  plane- 


parallel  atmospheric  layer  with  no  scattering  may  be  written 
I(x;u)  «  I(x1;u)e’(Tl'T)/li  +  /  Bv(T(t'  ))e-(x’  't)/y  ^  . 


The  cloud  is  assumed  isothermal  and  the  source  function,  B.^(T),  may  be 
removed  from  the  integral  and  the  solution  to  Eq.  may  be  written 

I (t  ;jj )  =  I (t i  ) e  '  +  Bv(T)[l-e  1  ]. 
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Note  that  at  90°,  the  solution  becomes 
I(t;u=0)  =  BV(T). 

Some  additional  physical  insight  may  be  gained  from  this  solution.  If 
the  cloud  optical  depth  is  zero  then  the  second  term  in  the  solution 
goes  to  zero  and 

I(t;y)  =  I(x-|  ;u) . 

When  the  cloud  is  optically  thick  the  exponential  factor  approaches 
zero  and  the  non-scattering  solution  again  approaches  the  Planck  func¬ 
tion.  Multiple  scattering  may  be  deleted  from  the  CL  model,  also,  but 
the  problem's  complexity  is  not  appreciably  reduced.  The  equation  may 
be  written 

[sinecos^  .  ,<*£)  .  t(f£)  „  BjT)  . 

Assuming  the  same  expansion  as  before  for  the  intensity  we  arrive  at  a 
recursion  formula  similar  to  Eq.  (5-59).  Using  the  method  of  selecting 
the  zero  and  first  order  modes,  we  may  continue  with  the  same  method 
to  find  a  solution  for  Ig(p«<t>r)-  It  appears  that  the  major  difference 
in  the  solution  is  in  the  eigenvalues,  A.  With  scattering,  the  first 
order  approximation  had 

A  = 

Without  scattering  the  eigenvalue  for  the  first  order  approximation 


becomes 


a  constant. 


The  numerical  program  was  modified  to  remove  the  multiple  scatter¬ 
ing  term  and  Fig.  23  shows  graphically  the  reduction  in  the  fifth  order 
solution  of  weighted  average  intensity  from  the  cloud  alone  at  a  fre¬ 
quency  of  50.50  GHz  (zero  outer  boundary  conditions)  when  multiple 
scattering  is  removed.  It  appears  that  the  effect  of  multiple  scatter¬ 
ing  is  to  increase  the  total  brightness  to  near  an  optical  depth  of  10. 
At  values  of  optical  depth  larger  than  10  the  curves  nearly  parallel 
each  other  with  multiple  scattering  contributing  a  near  constant  addi¬ 
tional  373  to  the  cloud  brightness. 

We  discuss  the  results  of  the  cloud  alone,  e.g.,  the  boundary  con¬ 
ditions  are  assumed  zero  and  the  cloud  is  given  a  temperature.  Upon 
examining  the  character  of  the  solution  at  several  different  points 
around  the  cylinder  cross  section,  it  was  found  that  the  emergent 
angular  intensity  distributions  were  identical.  This  demonstrates,  as 
we  might  expect,  that  the  cloud  with  zero  outer  boundary  conditions 
emits  isotropically. 

Non-zero  variably  boundary  conditions,  then,  will  act  as  a  forcing 
function  on  top  of  a  constant  solution.  The  amount  of  radiation  emit¬ 
ted  by  the  cloud  is,  of  course,  a  function  of  the  cloud  composition. 

Figures  24  and  25  present  the  variation  of  brightness  temperature 
with  emergent  angle  for  the  PP  and  CL  models  for  two  representative 
optical  depths  in  the  microwave  region.  Recall  that  in  the  PP  (dis¬ 
crete  ordinate)  method,  the  solution  is  gotten  explicitly  at  discrete 
emergent  angles  (Gauss  nodes)  and  then  a  curve  is  fit  through  those 
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Figure  24.  50.50  GHz  cloud  brightness  temperature  versus  emergent 

angle:  zero  outer  boundary  conditions.  Tcioud=^01° 
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to  extrapolate  the  0°  and  90°  solutions.  This  method  tends  to  present 
a  smooth  curve.  The  solution  as  described  for  the  cylinder,  at  one 
point,  is  a  series  of  spherical  harmonic  functions  and  will  tend  to 
have  some  oscillation,  depending  on  the  order  of  approximation  used. 

Dave  and  Canosa  (1974)  presented  results  for  a  spherical  harmonics 
approximation  to  plane-parallel  atmospheric  models.  Their  figures 
show  a  pronounced  oscillation  for  a  fifth  order  (so-called  Pg)  approxi¬ 
mation  for  optical  depths  up  to  1 .  The  various  orders  of  solutions 
for  t=1.8  for  the  CL  model  have  no  oscillation  and  we  begin  to  see 
some  oscillation  for  the  r=12.  curves.  This  slow  onset  of  oscillatory 
behavior  may  be  attributed  to  the  form  of  the  two  dimensional  solution. 
For  example,  please  refer  to  Eq.  5-82,  the  3rd  order  approximation 
written  out  explicitly.  The  solution  is  expressed  in  terms  of  various 
combinations  of  spherical  harmonics  (Eq.  5-57),  all  of  which  did  not 
vanish  when  the  emergent  angle  is  equal  to  90°.  Therefore,  the  approxi¬ 
mation  renders  a  value  at  90°  which  is  other  than  zero.  Also,  the  com¬ 
bination  of  the  functions  tends  to  damp  the  oscillation  experienced  by 
Dave  and  Canosa  for  low  orders  of  approximation  for  plane  parallel 
media . 

Cloud  emissivity  for  50.50  GHz  is  presented  in  Fig.  25.  The  aver¬ 
age  emissivity  values,  defined  as  the  ratio  of  the  average  outward 
flux  denisty  to  ttB.^T  ),  are  presented  as  a  function  of  cloud  optical 
thickness  for  first,  third  and  fifth  order  approximations  to  the  cylin¬ 
der  cloud  model.  The  emissivity  curves  have  essentially  the  same  trends 
as  presented  in  the  infrared  case  study  (Fig. 10).  Note,  however,  that 
the  values  are  lower  in  the  microwave  region  than  in  the  infrared  area 
for  both  models  at  very  large  optical  thicknesses.  Also  the  fifth 
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order  CL  solution  Is  less  than  the  plane-parallel  solution,  lending 
more  credence  to  the  conjecture  that  the  more  finite  a  cloud  model  Is, 
the  lower  Its  emlsslvlty  will  be. 


CHAPTER  7 


ANALYSIS  OF  DMSP  SSM/T  DATA 

7.1  Cases  Selected  for  Study 

We  have  selected  two  days,  30  October,  1979  and  23  November, 
1979  during  which  significant  cloud  and  precipitation  events  are 
present  and  during  which  complete  data  sets  were  available  over 
the  continental  United  States,  The  data  included  synoptic  re¬ 
ports  and  charts,  radar  summaries,  radiosonde  observations,  DMSP 
infrared  satellite  pictures  and  SSM/T  microwave  sounder  data. 

The  SSM/T  instrument  is  on  board  a  sun-synchronous  polar  orbiting 
satellite. 

The  cases  selected  are  presented  in  Table  6  and  the  satellite 
pass  times  presented  in  Table  6  are  the  actual  observation  times 
of  the  SSM/T  instrument.  Note,  that,  with  one  exception,  the  scan 
times  are  between  0000  Zulu  and  0600  Zulu  on  the  respective  dates. 
Therefore  the  0000  Zulu  radiosonde  observations  were  deemed  most 
representative  and  thus  were  used  in  the  following  studies.  The 
DMSP  infrared  satellite  pictures  for  30  October  79  are  from  the 
same  spacecraft  as  the  SSM/T  data  and  therefore  are  sensed  at  the 
same  time.  The  infrared  pictures  for  the  23  November  79  case  study 
have  descending  node  crossing  times  of  0043  Zulu  for  the  eastern 
picture  and  0224  Zulu  for  the  western  picture.  Figures  27  ,  28  and 
29  present  the  satellite  pictures  used  in  the  case  studies.  A  major 
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Table  6-  Selected  Cases 


Station 

Name 

Satellite 
Pass  Time 

Latitude 

(°N) 

Longitude 

(°W) 

Case  Type 

Jackson 

Mississippi  (MS) 

0346  Z 

30  Oct  79 

32.19 

90.05 

Clear 

Centerville 

Alabama  (AL) 

0346  Z 

30  Oct  79 

32.54 

87.15 

Clear 

Little  Rock 

Arkansas  (AR) 

0346  Z 

30  Oct  79 

34.44 

92.14 

Clear 

Monett 

Missouri  (MO) 

0347  Z 

30  Oct  79 

36.53 

93.54 

Clear 

Salem 

Illinois  (IL) 

0347  Z 

30  Oct  79 

38.39 

38.58 

Clear 

Peoria 

Illinois  (IL) 

0348  Z 

30  Oct  79 

40.40 

89.41 

Clear 

Boise 

Idaho  (ID) 

0530  Z 

30  Oct  79 

43.34 

116.13 

Clear 

Desert  Rock 

Nevada  (NV) 

0529  Z 

30  Oct  79 

36.6 

116.1 

Clear 

Greensboro 

North  Carolina  (NC) 

0251  Z 

23  Nov  79 

36.03 

79.57 

Clear 

Chihuahua 

Mexico  (MEX) 

0430  Z 

23  Nov  79 

28.42 

106.04 

Clear 

Del  Rio 

Texas  (TX) 

0430  Z 

23  Nov  79 

29.22 

100.55 

Clear 

Dodge  City 

Kansas  (KS) 

0432  Z 

23  Nov  79 

37.46 

99.58 

Clear 

Great  Falls 

Montana  (MT) 

0435  Z 

23  Nov  79 

47.29 

111. 22 

Clear 

Glasgow 

Montana  (MT) 

0435  Z 

23  Nov  79 

48.13 

106.37 

Clear 

Bismark 

North  Dakota  (ND) 

0350  Z 

30  Oct  79 

46.46 

100.45 

Cloudy 

Salem 

Oregon  (OR) 

0531  Z 

30  Oct  79 

44.55 

123.01 

Cloudy 

Medford 

Oregon  (OR) 

0531  Z 

30  Oct  79 

42.22 

122.52 

Cloudy 

Huntington 

W.  Virginia  (WV) 

0252  Z 

23  Nov  79 

38.22 

82.33 

Cloudy 
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Table  6.  Selected  Cases  (Continued) 


Station 

Name 

Satellite 
Pass  Time 

Latitude 

(°N) 

Longitude 

(°W) 

Case  Type 

Green  Bay 

Wisconsin  (WI) 

0253  Z 

23  Nov  79 

44.29 

88.08 

Cloudy 

Monterrey 

Mexico  (MEX ) 

0429  Z 

23  Nov  79 

25.52 

100.12 

Cloudy 

Victoria 

Texas  (TX) 

0429  Z 

23  Nov  79 

28.51 

96.55 

Cloudy 

Rapid  City 

S.  Dakota  (SD) 

0434  Z 

23  Nov  79 

44.03 

103.04 

Cl oudy 

Medford 

Oregon  (OR) 

0614  Z 

23  Nov  79 

42.44 

122.52 

Cl oudy 

Dodge  City 

Kansas  (KS) 

0348  Z 

30  Oct  79 

37.46 

99.58 

Precipitating 

Omaha 

Nebraska  (NE) 

0349  Z 

30  Oct  79 

41.22 

96.01 

Precipitating 

Huron 

S.  Dakota  (SD) 

0350  Z 

30  Oct  79 

44.23 

98.13 

Precipitating 

Dayton 

Ohio  (OH) 

0252  Z 

23  Nov  79 

39.52 

84.07 

Precipitating 

Pittsburgh 
Pennsylvania  (PA) 

0252  Z 

23  Nov  79 

40.32 

80.14 

Precipitating 

Buffalo 

New  York  (NY) 

0252  Z 

23  Nov  79 

42.56 

78.44 

Precipitating 

Spokane 

Washington 

0435  Z 

23  Nov  79 

47.38 

117.32 

Precipitating 

Figure  27.  30  Oct  79  DMSP  satellite  picture  mosaic 


Figure  28.  23  Nov  79  DMSP  satellite  -  eastern  scan 
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criterion  for  case  selection  was  to  have  a  radiosonde  station  with¬ 
in  the  microwave  sounder  field  of  view.  The  satellite  pictures 
have  been  superimposed  with  the  SSM/T  sounder  "footprints"  and  the 
relative  location  of  the  radiosonde  stations  within  them. 

7.1.1  Brief  Synoptic  Discussion  -  30  October,  1979 

A  surface  low  pressure  center  (Fig.  30)  lies  off  the  east  coast 
of  the  United  States  near  44°N.  62°W.  A  cold  front  extends  south¬ 
ward  and  then  west  to  southern  North  Carolina,  becoming  a  warm  front 
in  northern  Georgia  to  southern  Arkansas  and  into  Oklahoma.  A  dis¬ 
sipating  low  aloft  is  portrayed  on  the  30/0000  Zulu  500  mb  chart 
(Fig.  31)  as  a  weakness  through  the  Ohio  valley. 

The  dominant  feature  on  the  500  mb  map  is  a  storm  system  that 
is  deepening  and  moving  into  northeast  Arizona.  A  surface  low  pres¬ 
sure  center  is  drawn  in  southeast  New  Mexico.  Extending  northward  from 
this  surface  low  is  an  inverted  trough  extending  north-northeast  to 
central  Minnesota.  The  500  mb  chart  also  shows  this  trough  which  ex¬ 
tends  to  a  low  center  near  64°N,  107°W.  Cloud  patterns  associated 
with  the  southern  storm  are  extensive  (Fig.  27  )  and  case  studies  were 
chosen  in  and  around  this  system.  Mote  that  the  satellite  footprints 
and  selected  station  locations  are  presented  on  the  satellite  picture 
(Fig.  27).  The  details  of  each  individual  case  will  be  presented  in 
Section  7.2.  A  trough  at  500  mb  just  off  the  west  coast  of  the  U.S. 
supports  a  surface  low  drawn  near  49°N,  136°W.  An  occlusion  extends 
south-southeast  to  a  triple  point  near  43°N,  134°W.  This  system  is 
moving  rapidly  and  precipitation  begins  at  stations  on  the  Washington 


Figure  30.  30  Oct  79,  OOOOZ,  surface  analysis. 
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and  Oregon  coasts  ny  0600  Zulu.  The  cloud  pattern  associated  with 
this  system  appears  to  be  well  organized  with  the  cirrus  shield 
already  on  shore  at  the  satellite  pass  time  (Fig.  27). 

7.1.2  Brief  Synoptic  Qiscussion  -  23  November,  1979 

A  dominant  500  mb  low  pressure  center  (Fig.  32)  is  located  in 
central  Minnesota  at  23/0000  Zulu.  A  deep  500  mb  trough  extends  south- 
southwest  to  the  Texas  panhandle  and  continues  southwest  into  northwest 
Mexico.  On  the  surface  (Fig.  33)  there  is  a  double  low  system  with 
one  low  center  in  south  central  Minnesota  with  an  occluded  front  to  a 
second  low  drawn  on  the  triple  point  in  central  lower  Michigan.  A 
warm  front  continues  eastward  to  extreme  southern  Maine.  A  cold  front 
extends  south  through  Indiana  and  into  the  Gulf  of  Mexico  near  New 
Orleans,  Louisiana.  Considerable  cloudiness  accompanies  this  storm 
system  (Fig.  28),  extending  from  the  Canadian  border  across  Mexico  and 
into  the  eastern  Pacific.  Another  system  is  just  off  the  west  coast 
of  the  U.S.  at  500  mb  (Fig.  32).  The  surface  map  (Fig.  33)  shows  a 
cold  front  just  on  the  west  coast  from  northern  California  through 
western  Oregon  and  Washington.  Extensive  cloudiness  also  accompanies 
this  storm  (Fig.  29)  with  cloud  bands  extending  along  the  west  coast  of 
the  U.S. 

7.2  Comparison  of  Computed  and  Observed 

Brightness  Temperature 

The  case  studies  were  divided  into  two  groups  -  clear  and  cloudy/ 
precipitating.  The  clear  atmosphere  upwelling  brightness  temperatures 
were  computed  by  the  olane-parallel  method.  The  cloudy/precipitating 
cases  were  analyzed  using  the  plane-parallel  and  cylindrical  models. 


Figure  33.  23  Nov  79,  0000Z,  surface  analysis. 
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All  results  attributed  to  the  cylinder  model  in  the  following  sections 
are  fifth  order  solutions. 


7.2.1  Clear  Atmosphere  Cases 

A  major  problem  in  microwave  remote  sensing  over  land  is  determin¬ 
ing  the  surface  emissivity.  Gloersen,  et.  al .  (1972)  have  presented 
microwave  emissivity  values  ranging  from  approximately  .30  over  wet 
bare  soil  to  .98  over  a  uniformly  vegetated  surface.  This  background 
variability  coupled  with  atmospheric  uncertainties  (e.g.,  clouds  and 
precipitation)  considerably  complicate  the  remote  sensing  effort.  In 
the  clear  (non-scattering)  atmosohere  the  radiative  transfer  problem 
is  somewhat  simplified  and  the  surface  emissivity  may  be  estimated  by 
an  objective  scheme.  Two  such  schemes  were  investigated  for  use  in  the 
clear  atmosphere. 

The  first  method  involoves  solving  the  radiative  transfer  eguation 
for  the  surface  emissivity  using  the  window  frequency  of  50.50  GHz. 

This  frequency's  weighting  function  peaks  at  or  near  the  surface  and 
we  may  construe  that  the  bulk  of  the  energy  received  by  an  observer  at 
the  top  of  the  atmosphere  at  this  frequency  had  its  origin  at  the  sur¬ 
face.  This  basic  equation  may  be  written  as 


T  T 

S  Vj 


(0,»)  +  /  T(z)  dT  (z,») 
z=0  Jl 


:i-c  )  T  (Qr)  ;°  T(z)  dT  (0 ,z) 

Z=,o  vi 


(7.1) 


where  Tg  is  the  brightness  temperature  at  j i ( 50 .50  GHz),  o  is  the 

emissivity,  T  is  the  surface  temperature,  T 

s  J  j 


is  the  50.50  GHz 
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transmission  function,  and  T(z)  os  the  temperature  at  height  z.  For 
simplicity  we  may  define  the  following  terms: 


A-j=  Tvi  (0>°°) 

CO 

B, =  /  T(z)  dT  (z,«) 

1  z=0  Vl 

C, «  /°  T(z)  dT  (O.z) 

* 


(7.2) 


Substituting  Eq.  (7.2)  into  Eq. 

may  be  written 

c  =TB(v!)-  Br  A1C1 
- 

TsAl  -  A,C, 


(7.1)  and  solving  for  s,  the  equation 


(7.3) 


For  the  clear  column  case  studies,  the  atmosphere  radiosonde  pro 
files  may  be  used  to  generate  the  transmission  function,  the  surface 
temperature  averaged  over  the  satellite  "footprint"  may  be  estimated 
from  the  synoptic  reports,  and  the  observed  values  of  brightness 
temperature  may  be  used  for  Tg(vl).  Note  Eq.  (7.3)  may  be  solved  for 
the  50.50  GHz  surface  emissivity.  It  has  been  generally  found  that 
the  emissivity  does  not  change  in  value  between  50.50  and  59.40  GHz 
(Gloersen,  et.  al . ,  1972).  Thus,  the  50.50  GHz  value  may  be  employed 
in  upwelling  brightness  temperature  calculations  for  the  remaining 
channels.  Since  Eq.  (7.1)  is  applicable  to  all  frequencies,  we  may 
replace  v  (50.50  GHz)  in  Eq.  (7.1)  by  v;(53.20  GHz),  the  wing  channel 
whose  weighting  function  peak  was  nearest  the  surface.  The  equation 
for  50.50  GHz  was  then  divided  by  that  of  53.20  GHz  to  give 


109 


eiTsAr  Bi  "  Aici  +  eiAici 

£2TsA2=  TB2  '  B2  "  A2C2  +  £2A2C2 


(7.4) 


where  the  subscript  1  refers  to  50.50  GHz  and  2  refers  to  53.20  GHz. 
Once  again  the  surface  emissivity  is  assumed  invariant  for  the  fre¬ 
quency  range  considered  and  Eq.  (7.4)  may  be  solved  for  the  surface 
emissivity  independent  of  the  surface  temperature  as 


e  = 


fB1  "  B1  "  A1C1  ‘  ^2^TB2  ‘  B2  ’  A2C2) 
A1C2  ’  C1 


(7.5) 


This  emissivity  value  may  now  be  used  to  solve  Eq.  (7.1)  for  the  sur¬ 
face  temperature.  Surface  temperature  and  emissivity  values  generated 
by  methods  one  and  two  are  presented  in  Table  7.  The  emissivities 
calculated  by  method  1  seem  to  be  more  in  accordance  with  published 
values  (Gloersen,  et.  al . ,  1972)  for  the  different  types  of  ground 
surfaces  than  those  of  method  2.  Therefore,  the  emissivities  of  method 
1  were  used  for  the  clear  column  calculations  of  the  remaining  SSM/T 
frequencies.  These  calculations  are  presented  in  Fig.  34.  Note  that 
a  50.50  GHz  diagram  is  not  presented  since  the  calculated  and  ob¬ 
served  values  will  match  exactly  for  this  frequency  when  using  method 
1.  At  53.20  GHz  there  is  generally  good  agreement  between  the  cal¬ 
culated  and  observed  brightness  temperatures,  with  most  of  the  calcu¬ 
lations  slightly  underestimating  the  observed  values.  53.45  GHz  and 
54. 9C  GHz  channels  overestimate  the  observations  -  54.90  GHz  by  a 
larger  amount  in  the  average  than  54.35  GHz.  Referring  to  Fig.  3, 
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Table  7-  Generated  Emi ssivities  and  Surface  Temperatures 


Station 

Method  1 
e  Ts(°K) 

Method  2 
£  TS(°K) 

Jackson,  MS 

.954 

289.8 

.890 

303.4 

Centreville,  AL 

.949 

289.8 

.895 

301.1 

Little  Rock,  AK 

.958 

287.1 

.901 

299.4 

Monett,  MO 

.971 

287.1 

.902 

301  .9 

Salem,  I L 

.969 

282.1 

.873 

302.5 

Peoria,  IL 

.940 

285.9 

.894 

295.9 

Dodge  City,  KS 

.933 

273.2 

.396 

281.3 

Del  Rio,  TX 

.924 

284.6 

.905 

288.9 

Great  Falls,  MT 

.908 

275.9 

.933 

270.6 

Glasgow,  MT 

.902 

281 .4 

.927 

275.7 

Greensboro,  NC 

.935 

290.9 

.373 

305.2 

Chihuahua,  MEX 

.912 

280.6 

.948 

272.4 

Desert  Rock,  NV 

.900 

281.5 

.814 

307.5 

Boise,  ID 

.906 

273.2 

.921 

270.2 
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the  weighting  function  diagram,  it  is  seen  that  all  three  of  these 
frequencies'  weighting  functions  peak  well  below  16  km.  The  temp¬ 
erature  profile  is  known  to  be  at  least  16  km  in  all  of  the  clear 
cases.  In  theory,  then,  the  calculated  values  should  match  the 
observations  to  within  the  accuracy  of  the  radiosonde  observation. 

The  brightness  temperature  at  53.20  GHz  frequency  is  signifi¬ 
cantly  affected  by  the  surface.  Since  the  height  of  the  standard 
meteorological  instrument  shelter  is  6  feet,  tnis  should  cause  the 
actual  surface  temperature  to  be  underestimated,  This  fact  may  pro¬ 
vide  a  possible  answer  for  53.20  GHz  calculation  generally  being  lower 
than  the  observed  values. 

The  brightness  temperatures  at  54.35  GHz  and  54.90  GHz  are  not 
appreciably  affected  by  the  surface,  since  the  transmission  from  the 
Earth's  surface  to  the  top  of  the  atmosphere  is  approximately  .023 
and  .003  respectively.  Since  the  temperature  profile  is  known  to  a 
reasonable  accuracy  at  these  heights,  it  would  appear,  based  on  the 
studies,  that  there  may  be  some  systematic  error  involved  in  cal¬ 
culating  the  brightness  temperature.  This  systematic  overestimation 
in  remote  sensing  channels  has  been  noted  by  Liou  and  Yeh  (personal 
communication)  in  the  infrared  region.  Since  the  transmission  func¬ 
tions  are  calculated  in  a  different  manner  for  the  infrared  region 
than  in  the  microwave  region,  and  the  error  is  still  apparent,  it 
seems  that  the  possibility  of  faulty  transmission  functions  is  re¬ 
duced.  A  more  likely  reason  may  be  that  the  source  function  for  the 
clear  column  radiation  calculations  is  somehow  being  systematically 
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overestimated.  This  conjecture  was  first  mentioned  by  McClatchey 
et.  al.  (1976). 

Radiosonde  observations  were  not  available  above  16  km  for  the 
clear  column  case  studies  with  the  exception  of  Chihuahua,  Mexico  on 
23  November  1979,  where  the  observation  was  24  km.  Where  the  ob¬ 
servation  was  not  available,  the  Spring/Fall  climatological  profile 
was  substituted.  Referring  to  Fig.  3,  it  is  seen  that  the  weighting 
functions  of  frequencies  of  58.82  GHz,  59.40  GHz  and  58.40  GHz  peak 
above  16  km.  The  58.82  GHz  weighting  function  peaks  highest  in  the 
atmosphere  with  only  negligible  energy  emerging  from  the  top  of  the 
atmosphere  originating  below  18  km.  As  a  result,  the  calculated  values 
were  nearly  the  same  for  all  of  the  cases  since  they  were  all  based  on 
the  same  climatological  sounding.  58.82  GHz  and  59.40  GHz  weighting 
functions  peak  lower  in  the  atmosphere,  at  17  and  20  km,  respectively, 
and  should  be  influenced  somewhat  by  the  available  soundings.  It  is 
evident  that  the  plots  are  approaching  the  results  portrayed  for 
53.40  GHz  as  the  peak  of  the  weighting  function  rises  higher  in  the 
atmosphere  (e.g.,  the  calculated  values  are  more  and  more  representa¬ 
tive  of  the  climatological  sounding).  The  exception  to  this  is  the 
value  presented  for  Chihuahua,  Mexico,  whose  sounding  was  available  to 
25  km.  The  results  for  this  sounding  are  underlined  (Fig.  34)  on  the 
plot  for  58.82  GHz  and  59.40  GHz  and  also  demonstrates  that  the  calcu¬ 
lated  value  is  larger  than  the  observed  brightness  temperature,  repeat¬ 
ing  the  trend  established  by  the  comparisons  at  54.35  and  54.90  GHz. 
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7.2.2  Cloudy/Precipitating  Cases 

As  stated  in  Chapter  3,  the  Marshall -Palmer  drop  size  dis¬ 
tribution  was  used  to  simulate  precipitating  clouds.  Non-precip¬ 
itating  clouds  were  modeled  using  a  drop  size  distribution  established 
by  experimentation  and  presented  by  Mason  (1971),  as  measured  by  Diem 
(1948)  for  cumuliform  clouds. 

Figure  35  presents  how  the  cylinder  model  may  be  applied  to  the 
atmosphere  within  the  satellite  field  of  view.  This  sort  of  config¬ 
uration  suggests  a  cold  frontal  band  or  squall  line  in  the  real  at¬ 
mosphere.  The  SSM/T  instrument  chosen  for  case  study  investigations 
has  a  footprint  on  the  Earth's  surface  that  ranges  from  a  near  circle 
of  radius  of  87  km  at  nadir  to  an  ellipse  with  a  semi-major  axis  of 
152  km  at  maximum  scan  angle  (see  Fig.  2). 

In  order  to  realistically  evaluate  the  upwelling  brightness  temp¬ 
erature  in  cloudy  atmospheres  within  the  SSM/T  field  of  view,  modifi¬ 
cation  of  the  cylinder  model  is  required  so  that  cloud  and/or  pre¬ 
cipitation  could  be  placed  below  the  tropopause.  The  cylinder  was 
placed  so  that  a  portion  of  it  was  below  the  surface.  The  portion 
remaining  above  the  surface  was  dictated  by  the  height  of  the  cloud 
top  and  the  areal  extent  of  the  cloud  formation.  For  outer  boundary 
conditions,  incoming  intensity  streams  on  the  portion  of  the  cylinder 
beneath  the  surface  were  chosen  to  be  eT$ .  Since  the  size  of  the 
cylinder  will  be  increased,  the  amount  of  scattering  material  in  the 
cylinder  will  also  be  increased.  This,  in  effect,  would  cause  an  un¬ 
realistically  large  optical  depth,  due  to  scattering,  in  the  portion 
of  the  cylinder  that  represented  the  cloud.  For  this  reason,  adjust- 
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ment  was  made  such  that  the  optical  depth  due  to  scattering  was  re¬ 
duced  by  the  ratio  of  the  area  of  the  cylinder  cross  section  repre¬ 
senting  the  cloud  to  the  total  cross  sectional  area.  In  practice, 
this  ratio  was  multiplied  by  the  scattering  coefficient.  The  cloud 
temperature  was  chosen  as  the  midpoint  value  along  the  radius  vector 
at  nadir-zenith.  Finally,  since  only  a  portion  of  the  field  of  view 
may  be  covered  by  the  cloud,  the  cylinder  cloud  solution  is  averaged 
with  a  clear  column  solution  based  on  the  temperature  profile  at  the 
radiosonde  station  in  the  field  of  view.  The  solutions  are  weighted 
according  to  the  precentages  of  the  field  of  view  covered  by  cloud. 

Precipitating  and  cloudy  cases  are  presented  individually.  Be¬ 
cause  the  methods  used  for  deriving  surface  emissivity  in  the  clear 
atmosphere  could  not  be  applied  in  these  cases,  the  was  fixed  at 
.97,  a  compromise  value  based  on  published  literature  (Gloersen,  et. 
al.,  1972).  Note  that  the  precipitating  cases  are  in  areas  that  are 
heavily  vegetated,  and  as  a  result,  the  surface  emissivity  should  not 
degrade  significantly  during  precipitation.  Estimations  of  the  thick¬ 
ness  and  cloud  amount  within  the  field  of  view  were  made  using  the 
available  radiosonde  data  and  DMSP  infrared  satellite  photographs. 

The  plane  parallel  cloud  used  for  comparison  with  the  cylinder  cloud 
was  composed  of  three  isothermal  layers.  Because  of  the  plane  parallel 
model's  geometry,  the  plane  parallel  cloud  filled  the  entire  field  of 
view. 

The  information  gathered  for  input  to  the  cloudy  cases  is  pre¬ 
sented  as  follows: 

1.  Bismark,  N.O.  -  30  Oct  79.  The  satellite  picture  indicates  a 


117 


cirrus  shield  over  low  cloudiness.  Radar  charts  show  no  precipitation 
in  the  field  of  view.  Station  reports  and  the  National  Weather  Service 
nephanalysis  show  low  cloudiness  within  the  field  of  view  covering  an 
estimated  60%  of  the  scan  area  with  cloud  bases  at  .9  Km  and  tops  at 
2.75  Km.  All  bases  and  tops  given  are  above  ground  level  values.  The 
scan  angle  is  24°. 

2.  Medford,  OR.  -  30  Oct  79.  Satellite  picture  indicates  cloudi¬ 
ness  in  the  northwest  portion  of  the  scan  spot.  Synoptic  data  indicate 
this  is  a  band  of  alto-cumulus  with  bases  at  4.5  Km  and  top  at  6.0  Km. 
Cloud  coverage  is  estimated  30%  and  the  scan  angle  is  12°. 

3.  Salem,  OR.  -  30  Oct  79.  Satellite  picture  shews  overcast  cirrus. 
Synoptic  reports  and  nephanalysis  indicate  a  layer  of  alto-stratus  below 
the  cirrus  between  3.5  and  4.5  Km.  The  alto-stratus  is  estimated  to 
fill  40%  of  the  field  of  view.  The  scan  angle  is  12°. 

4.  Rapid  City,  S.D.  -  23  Nov  79.  Satellite  picture  shows  an  area 
pf  alto-cumulus  with  little  or  no  organization.  Synoptic  reports  tend 
to  confirm  this.  Bases  were  estimated  at  2  Km  and  tops  at  3  Km.  Cloud 
coverage  was  assumed  to  be  50%.  The  scan  angle  is  12°. 

5.  Huntington,  W.V.  -  23  Nov  79.  Sa  3 11 i te  picture  indicates 
cirrus  over  some  low  cloudiness.  Synoptic  reports  show  the  lower  cloud¬ 
iness  to  be  cumulus  and  cumulus  congestus.  Cloud  bases  are  assumed  to 
be  1.25  Km  amd  tops  at  4.0  Km  with  areal  coveraae  set  at  30%.  The  scan 
angle  is  12°. 

6.  Green  Bay,  WI.  -  23  Nov  79.  A  thin  cirrus  band  lies  over  a 
low  stratus  deck  on  the  backside  of  the  filling  surface  low  center. 
Synoptic  reports  indicate  the  cloud  bases  at  .1  Km  and  tops  at  1.6  'Km. 
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Coverage  is  estimated  at  90%  and  the  scan  angle  is  36°. 

7.  Monterrey,  MX  -  23  Nov.  79.  Satellite  picture  shows  a  large 
cloud  band  over  the  station  which  fills  most  of  the  scan  spot.  Synop¬ 
tic  reports  indicate  this  band  is  a  layer  of  alto-cumulus  with  ci rrus 
above.  Cloud  bases  were  assumed  to  be  2  km  with  tops  at  3  km.  The 
scan  is  at  0°  and  the  area  is  80%  cloud  covered. 

8.  Victoria,  TX  -  23  Nov.  79.  This  is  the  same  cloud  band  that 
is  over  Monterrey,  Mexico,  but  synoptic  reports  indicate  the  layer  of 
alto-cumulus  is  thicker,  with  bases  at  2  km  and  tops  at  4.5  km.  The 
scan  angle  is  24°.  Eighty  percent  cloud  coverage  was  found. 

Scatter  diagrams  are  presented  for  the  cloudy  cases  in  Figures 
36  and  37.  Note  that  both  models  give  brightness  temperature  values 
that  are  in  the  same  range.  Also,  calculations  from  both  models  con¬ 
sistently  give  values  larger  than  the  observations,  repeating  the 
trend  denoted  in  the  clear  column  study. 

The  precipitating  cases  are  presented  as  follows: 

1.  Dodge  City,  KS  -  30  Oct.  79.  Satellite  picture  indicates  a 
good  squall  line  in  the  field  of  view.  Since  the  scan  angle  is  36°, 
the  sensor  may  have  picked  up  some  signature  from  the  leading  cloud 
band  on  the  eastern  edge  of  the  satellite  footprint.  The  radiosonde 
station  appears  to  be  in  the  clear  at  the  scan  time  and  was  definitely 
in  the  clear  at  radiosonde  time.  The  cloud  band  thickness  was  esti¬ 
mated  to  be  from  the  surface  to  12  km  and  covered  30%  of  the  scan  area. 
Rainfall  rate  was  estimated  at  5  mm/hr. 

2.  Omaha,  NE  -  30  Oct.  79.  The  satellite  picture  indicates 
cloudiness  to  tne  northwest  and  southwest  of  Omaha.  This  situation 
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CALCULATED  BRIGHTNESS  TEMPERATURE  (*K) 


Figure  36.  53.  50  and  53.23  GHz  calculated  versus  observed  cloudy 

Drightness  temperatures . 


OBSERVED  BRIGHTNESS  TEMPERATURE  (°K) 


CALCULATED  BRIGHTNESS  TEMPERATURE  (®K) 


igure  37.  54.35,  54.90,  58.40,  58.82  and  59.40  GHz  calculated 

versus  observed  cloudy  case  brightness  temperature. 

:= . 37 ,  *  plane-parallel,  "x"  =  5th  order  cylinder 
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was  modeled  by  a  cylinder  covering  30%  of  the  scan  footprint.  Cloud 
base  was  assumed  to  be  the  surface  and  cloud  top  was  at  3  km.  The 
scan  angle  was  24°  and  the  rainfall  rate  was  assumed  to  be  5  mm /hr. 

3.  Huron,  SD  -  30  Oct.  79.  The  satellite  picture  shows  an 
extensive  cloud  area  in  the  scan.  Radar  observations  are  revealing, 
showing  the  eastern  portion  of  the  scan  to  have  no  precipitation.  The 
area  covered  by  precipitation  was  assumed  to  be  70%.  Cloud  tops  and 
bases  were  from  surface  to  5  km.  The  scan  angle  was  12°  and  the  rain¬ 
fall  rate  estimated  to  be  1  mm/hr. 

4.  Spokane,  WA  -  23  Nov.  79.  The  satellite  picture  indicates 
two  cloud  bands  in  the  field  of  view.  The  leading  cloud  band,  from 
synoptic  reports,  is  believed  to  be  ci rrus  with  some  alto-cumulus  be¬ 
low.  The  trailing  cloud  band  has  precipitation.  The  surface  tempera¬ 
ture  is  34°F  and  it  is  snowing  at  Spokane.  The  cloud  modeled  had  a 
base  at  the  surface  and  the  top  was  assumed  to  be  5  km.  The  area 
covered  with  precipitating  cloud  was  50%  of  the  total.  Rainfall  rate 
was  1  mm/hr  and  the  scan  angle  was  36°. 

5.  Dayton,  OH  -  23  Nov.  79.  The  satellite  picture  shows  cloud 
covering  nearly  all  of  the  scan.  Synoptic  data  reveals  that  the 
eastern  portion  of  the  cloudiness  may  be  cirrus  only.  The  raining 
cloud  thickness  was  assumed  to  be  from  the  surface  to  3.5  km  and  to 
cover  o0%  of  the  scan.  The  rainfall  rate  was  1  mm/hr  and  the  scan 
angle  was  24°. 

6.  Pittsburgh,  PA  -  23  Nov.  79.  Satellite  pictures  reveal  that 
clouds  completely  fill  the  field  of  view.  Radar  reports  indicate  that 
the  extreme  eastern  and  western  edges  of  the  scan  may  be  cirrus  only. 
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Therefore,  the  raining  cloud  was  assumed  to  cover  80%  of  the  scan.  The 
cloud  base  was  at  the  surface  and  tops  were  assumed  to  be  at  4.0  km. 
Rainfall  rate  was  1  mm/hr  and  the  scan  angle  was  0°. 

7.  Buffalo,  NY  -  23  Nov.  79.  Satellite  picture  shows  several 
showers  in  the  scan  but  no  clear  organization.  Radar  reports  confirm 
this.  Raining  cloud  coverage  was  assumed  to  be  40%  with  cloud  tops 
at  4  km.  The  scan  angle  was  12°  and  the  rainfall  rate  was  1  mm/hr. 

Scatter  diagrams  are  presented  for  the  precipitating  cases  in 
Figures  38  and  39.  The  diagrams  for  50.50  GHz  and  53.20  GHz  show  some 
interesting  results.  Generally  the  cylinder  model  predicts  brightness 
temperature  values  that  are  slightly  higher  than  observed.  The  plane- 
parallel  results  are  all  lower  than  the  observed  values  for  these  fre¬ 
quencies.  In  two  cases,  the  plane-parallel  results  were  considerably 
less  than  the  observed  while  the  cylinder  results  were  quite  good. 

This  characterise  c  continues  to  the  frequency  with  the  next  higher 
weighting  function  peak,  54.  35  GHz.  These  cases  correspond  to  the  only 
two  case  studies  that  had  rainfall  rates  greater  than  1  mm/hr.  The 
plane-parallel  model  failed  the  worst  at  Dodge  City,  KS.  The  other 
poor  case  was  Omaha,  ME.  Both  cases  had  5  mm/hr  rainfall  rates.  Ac¬ 
tual  rainfall  rates  were  determined  by  estimating  the  average  rainfall 
in  a  scan  from  the  observed  rainfall  amounts  and  then  dividing  that 
value  by  the  time  the  radar  echoes  remained  over  the  area. 

Some  additional  insight  may  be  gained  from  these  studies.  The 
Spokane,  WA  case  supports  an  argument  put  forth  by  Savage  (1976)  that 
snowflakes  (as  opposed  to  smaller  ice  crystals  found  in  cirri  form 
clouds  which  have  negligible  effect  on  these  microwave  frequencies) 
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Figure  38.  50.50  and  53.20  GHz  calculated  versus  observed  precipi¬ 
tating  case  brightness  temperatures.  e=.97,  = 

plane-parallel,  "x"  =  5th  order  cylinder. 
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CALCULATED  BRIGHTNESS  TEMPERATURE  (*K) 


Figure  39.  54.35,  54.90,  58.40,  58.82  and  59.40  GHz  calculated  versus 

observed  precipitating  case  brightness  temperatures. 
e=.97,  =  plane-parallel,  "x"  =  5th  order  cylinder. 
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have  a  pronounced  effect  on  microwave  radiative  transfer.  Also,  the 
question  arises  as  to  why  the  plane-parallel  model  failed  so  badly  at 
the  higher  rainfall  rates.  Inspection  of  the  cylinder  cloud  solution 
(not  weighted  with  the  clear  coluim  for  area  coverage  consideration) 
shows  values  of  243.75  and  252.99°  K  for  Dodge  City,  KS  and  Omaha,  NE, 
respecti vely,  as  opposed  to  the  218.29  and  234.3°  K.  values  for  the 
plane-parallel  model.  We  note  that  in  these  two  cases,  the  areal 
coverages  by  precipitation  were  only  30%.  Thus,  it  seems  that  a 
finite  cylindrical  model  gives  a  much  better  approximation  than  the 
plane-parallel  model  in  the  computations  of  the  upwelling  brightness 
temperatures  in  these  cases.  Tables  8  and  9  present  the  numerical 
results  that  are  depicted  in  Figures  36  through  39. 

7.3  Temperature  Retrieval 

The  prime  objective  of  the  DMSP  SSM/T  microwave  sounders  has  been 
to  derive  accurate  temperature  profiles  in  all  weather  conditions  in 
conjunction  with  numerical  weather  prediction.  In  the  previous  com¬ 
parison  program  for  the  computed  and  observed  brightness  temperatures, 
temperature  profiles  for  selected  synoptic  cases  have  been  obtained 
from  the  available  radiosondes.  Thus,  it  seems  interesting  and  signi¬ 
ficant  to  examine  the  accuracy  and  validity  of  atmospheric  conditions 
over  land.  In  this  study,  we  use  a  statistical  method  for  temperature 
recovery  intended  for  operational  use  by  the  Air  Force  Global  Weather 
Central . 

In  the  statistical  approach,  it  is  generally  assumed  that  the  de¬ 
rivation  of  the  predicted  parameter  from  the  climatological  mean  may 
be  expressed  as  a  linear  combination  of  the  deviation  of  the  measured 


126 


Table  8.  Comparisons  between  observed  and  calculated 

brightness  temperatures  for  precipitation  cases. 


30  Oct 

79 

1 

2 

3 

Channels 

4 

5 

6 

7 

06S 

264.81 

253.93 

233.59 

224.04 

229.10 

214.74 

217.15 

72451 

CL 

265.35 

254.47 

238.89 

230.66 

226.39 

218.72 

219.95 

PP 

218.29 

219.30 

221.40 

222.18 

226.56 

218.73 

219.75 

08S 

268.02 

258.85 

238.08 

226.97 

225.56 

213.63 

215.38 

72553 

CL 

266.36 

253.89 

238.5 

229.53 

226.72 

217.85 

219.58 

PP 

234.3 

234.57 

231.57 

227.38 

226.56 

217.85 

219.47 

OBS 

266.12 

255.47 

236.93 

226.7 

224.68 

215.66 

216.16 

72654 

CL 

268.82 

261.67 

237.93 

228.08 

226.75 

218.18 

219.65 

PP 

259.27 

251.97 

237.39 

228.67 

226.56 

218.15 

219.55 

23  Nov 

79 

OBS 

257.55 

249.24 

230.83 

222.77 

224.26 

219.05 

219.49 

72785 

CL 

258.45 

250.81 

234.16 

224.45 

227.15 

218.15 

220.03 

PP 

253.72 

248.49 

236.17 

227.48 

226.56 

217.87 

219.68 

OBS 

272.26 

259.92 

236.66 

223.75 

223.85 

209.89 

211.86 

72429 

CL 

276.99 

265. 14 

240.16 

229.45 

226.99 

218.31 

219.90 

PP 

268.15 

258.85 

241.39 

231.15 

226.56 

218.29 

219.73 

OBS 

269.30 

260.89 

239.59 

226.19 

222.06 

210.18 

210.91 

72520 

CL 

271.63 

267.88 

241.25 

228.56 

226.56 

215.53 

218.92 

PP 

268.20 

258.12 

239.46 

228.37 

226.56 

215.53 

218.92 

OBS 

270.62 

260.71 

238.23 

225.80 

222.42 

209.84 

211.26 

72528 

CL 

274.57 

263.38 

239.64 

228.37 

226.65 

216.3 

219.18 

PP 

265.10 

256.46 

239.30 

228.80 

226.56 

216.22 

219.09 
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Table  9.  Comparisons  between  observed  and  calculated 
brightness  temperatures  for  cloudy  cases. 


30  Oct 

79 

1 

2 

3 

Channels 

4 

) 

5 

6 

7 

OBS 

264.79 

253.32 

234.78 

225.75 

224.99 

216.94 

217.19 

72764 

CL 

266.38 

254.0 

235.29 

227.26 

226.88 

218.89 

219.91 

PP 

265.92 

253.33 

236,37 

228.20 

226.56 

218.89 

219.75 

OBS 

262.78 

255.31 

237.  11 

225.69 

227.39 

214.89 

216.6 

72597 

CL 

266.84 

256.51 

238.68 

228.75 

226.64 

217.41 

219.47 

PP 

266.99 

256.35 

238.96 

229. 18 

226.56 

217.37 

219.42 

OBS 

260.65 

255.32 

236.84 

225.68 

225.95 

215.47 

217.13 

726  9  4 

CL 

270.94 

260.93 

237.79 

227.34 

226.67 

217.48 

219.60 

PP 

267.69 

255.41 

237.77 

227.95 

226.56 

217.42 

219.53 

23  Nov 

79 

OBS 

255.58 

251.59 

234.78 

226.29 

224.03 

216.58 

216.07 

72662 

CL 

264.94 

254.69 

236.17 

227.89 

226.65 

219.1 

219.87 

PP 

264.14 

252.25 

235.87 

228.10 

226.56 

219.11 

219.81 

OBS 

272.88 

261.60 

239.34 

225.72 

224.09 

209.68 

211.33 

72425 

CL 

279.22 

263. 19 

240.52 

229.05 

226.75 

215.92 

219.15 

PP 

278.85 

263.39 

241.54 

229.71 

226.56 

215.76 

218.96 

OBS 

254.67 

251.69 

233.50 

224.03 

223.07 

213.95 

214.04 

72  6  45 

CL 

266.79 

254.86 

238.04 

229.37 

226.56 

218.04 

219.55 

PP 

266.71 

254.88 

238. 17 

229.41 

226.56 

218.04 

219.55 

OBS 

266.62 

258.86 

239.20 

227.09 

229.67 

210.43 

213.76 

76394 

CL 

274.99 

265.63 

240.94 

229.53 

226.56 

215.58 

218.8 

PP 

273.80 

259.59 

239.55 

229.12 

226.56 

215.58 

218.79 

OBS 

259.51 

257.86 

236.37 

224.81 

230.39 

210.93 

214.66 

72255 

CL 

268.23 

259.56 

240.49 

230.25 

226.99 

218.33 

219.92 

PP 

268.51 

258.37 

241.62 

231.79 

226.56 

218.22 

219.72 
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data.  Upon  finding  a  linear  operator  which  will  yield  a  minimum  mean 
square  deviation  of  the  predicted  temperature  profile  from  the  true 
temperature  profile  in  a  statistical  sense,  the  predicted  temperature 
profile  may  be  obtained.  The  linear  operator,  called  the  predictor 
matrix,  may  be  expressed  in  terms  of  a  covariance  matrix,  which  can  be 
constructed  experimentally  by  collecting  coincidences  of  radiances  de¬ 
rived  from  remote  sounders  with  temperature  values  obtained  from  di¬ 
rect  soundings.  In  the  microwave  region,  where  the  surface  emissivity 
varies  significantly  with  land  and  ocean,  it  is  desirable  to  remove 
its  effect  on  the  temperature  retrieval  program  so  that  the  predictor 
matrix  could  be  constructed  over  all  surfaces.  This  can  be  done  by 
utilizing  the  surface  channel  (Rigone  and  Strogryn,  1977). 

Preliminary  results  using  the  statistical  method  over  ocean  sur¬ 
faces  are  quite  good  (Captain  Joe  Gahlinger,  personal  communication) 
and  a  similar  method  is  being  research  by  the  National  Weather  Service. 
Temperature  profile  retrievals  over  land,  however,  continue  to  pose 
problems.  As  shown  by  Nipko  (1979)  using  theoretical  data,  the  re¬ 
trieved  temperature  profile  suffers  severe  degradation  in  the  presence 
of  precipitating  clouds.  Precipitating  areas  may  be  eliminated  quite 
easily  over  ocean  surfaces.  Recall  Fig.  12,  the  50.50  GHz  frequency 
over  a  simulated  ocean  surface  with  e=.51.  If  the  window  channel 
brightness  temperature  is  higher  than  a  certain  threshold  value,  there 
must  be  a  precipitating  cloud  within  the  field  of  view.  Over  land 
surfaces,  partly  because  of  the  variability  in  the  surface  emissivity, 
it  is  difficult  to  discern  which  clouds  are  precipitating  and  which 
are  not.  An  additional  point  that  should  be  made  is  that  generally 
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meteorologists  and  their  numerical  models  would  prefer  to  have  tempera¬ 
ture  profile  information  in  the  vicinity  of  active  weather  producing 
storms,  precisely  the  areas  that  appear  to  pose  difficulties.  With 
these  thoughts  in  mind,  the  statistical  method  was  applied  to  the  known 
clear,  cloudy  and  precipitating  cases  presented  earlier  in  this  chapter 
in  hope  of  shedding  further  light  on  the  temperature  profile  recovery 
problem  over  land. 

The  statistical  method  was  applied  to  the  cases  studied  and  the 
results  are  presented  in  Figs.  40  through  42  for  four  clear,  cloudy  and 
precipitating  cases.  The  model  yields  atmospheric  temperatures  for 
mandatory  pressure  levels  of  1000,  850,  700,  500,  400,  300,  250,  150, 
100,  70,  50,  30,  20,  and  10  mb,  thicknesses  between  the  mandatory 
levels,  and  the  temperature  and  pressure  of  the  tropopause.  The 
heights  of  the  mandatory  levels  are  computed  by  stacking  the  thick¬ 
nesses  onto  a  forecast  1000  mb  height  field  retrieved  from  the  AFGWC 
computer  system.  In  these  studies  the  0000  Zulu  height  field  on  the 
appropriate  day  was  used  for  stacking  purposes. 

Figure  40  shows  the  retrieved  temperature  profiles  for  the  four 
clear  cases;  two  on  30  October  1979  and  two  on  23  November  1979. 

Except  for  the  Greensboro  case,  the  retrieved  temperature  profiles 
when  they  are  compared  with  nearby  radiosonde  data  appear  to  be  rea¬ 
sonably  good  in  view  of  the  statistical  method  used.  The  failure  in 
the  retrieval  for  the  Greensboro  case  seems  largely  due  to  the  fluctu¬ 
ated  temperature  profile  that  occurred  in  the  atmosphere.  Generally, 
we  found  that  the  statistical  method  is  working  properly  when  the 
actual  profile  is  smooth  and  when  no  inversion  is  present.  The  re¬ 
trieved  temperature  profiles  under  cloudy  conditions  are  illustrated 


Figure  40.  AFGWC  statistical  temperature  retrievals  -  clear  cases 
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Figure  41.  AFGWC  statistical  temperature  retrievals  -  cloudy  cases. 
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in  Fig.  41.  Apparently,  the  temperature  profile  retrieval  program 
using  the  microwave  frequencies  in  the  60  GHz  oxygen  band  is  affected 
insignificantly  by  non-precipitating  clouds.  Compared  with  the  tem¬ 
perature  obtained  from  radiosonaes,  the  retrieved  patterns  involving 
clouds  are  similar  to  those  under  clear  conditions.  Mote  that  in  each 
diagram,  the  percentage  of  cloud  cover  is  depicted.  Generally,  the 
retrieved  and  observed  (radiosonde)  temperatures  in  clear  and  cloudy 
conditions  are  within  about  5°  K.  In  the  final  figure  (Fig.  42)  we 
show  the  retrieved  temperature  profiles  under  precipitating  conditions. 
Again,  four  cases  are  presented  in  this  study.  Two  cases  are  selected 
from  30  October  1979;  both  have  a  5  mm/hr  rainfall  rate  with  30?£  cloud 
cover  in  the  field  of  view  of  the  SSM/T.  In  the  other  two  cases, 
selected  from  23  Movember  1979,  both  indicate  a  1  mm/hr  rainfall  rate 
but  with  cloud  covers  varying  from  50"  to  80®.  The  most  distinct 
feature  in  the  retrieved  temperature  profiles  using  the  statistical 
covariance  method  for  precipitating  cases  is  the  significant  and  con¬ 
sistent  deviation  from  the  radiosonde  data  in  the  lower  boundary  layer 
where  preci pi tati on  takes  place.  In  the  moderate  5  mm/hr  rainfall  rate 
cases,  the  differences  between  the  retrieved  and  radiosonde  temperature 
profiles  near  the  surface  are  as  large  as  10  -  15°  K.  It  should  be 
noted  that  precipitation  in  these  two  cases  cover  only  about  30"  within 
the  field  of  view  of  SSM/T.  As  for  the  cases  involving  1  mm/hr  rain¬ 
fall  rate,  about  5  -  10°  K  differences  near  the  ground  are  observed. 

'he  findings  for  these  precipitation  cases  using  the  real  SSM/T 
data  are  in  general  agreement  with  those  described  by  Mipko  (1979)  who 
oerformed  hypothetical  temperature  retrieval  exercises.  In  Mipko's 
analyses,  a  p lane-parallel  microwave  radiative  transfer  model  was  used 


134 


to  simulate  the  brightness  temperature  at  the  top  of  the  atmosphere. 

He  assumed  a  2  km  thick  cloud  with  rainfall  rates  ranging  from  1  mm/hr 
to  30  mm/hr  and  used  the  simulated  brightness  temperatures  from  seven 
SSM/T  channels  to  recover  the  temperature  profile  also  utilizing  the 
statistical  covariance  method.  The  surface  temperatures  are  shown  to 
suffer  increased  degradation  from  10°  to  20°  K  as  the  rainfall  rate 
increases  from  1  mm/hr  to  10  mm/hr.  The  hypothetical  analysis  that  he 
has  carried  out,  of  course,  assumes  that  the  cloud  covers  the  entire 
field  of  view. 

Although  the  current  study  employs  only  four  precipitating  cases 
in  the  analysis  and  may  not  be  conclusive  in  view  of  the  limited  sam¬ 
ple  used,  it  appears  that  the  effect  of  precipitation  on  the  tempera¬ 
ture  profile  retrieval  using  microwave  frequencies  is  substantial  and 
significant.  Of  course,  the  reliability  of  the  statistical  method  for 
the  temperature  profile  retrieval  in  clear  atmospheres  should  be  exam¬ 
ined  comprehensively  and  completely  utilizing  the  data  that  are  avail¬ 
able  in  different  localities  and  seasons.  In  addition,  in  order  to 
derive  the  temperature  field  over  the  global  space,  further  studies 
concerning  the  influence  of  precipitating  clouds  on  the  temperature 
retrieval  in  the  microwave  region  seem  also  warranted. 


CHAPTER  8 


CONCLUSION 

A  two-dimensional  cylindrical  cloud  model  was  constructed  for 
the  purpose  of  investigating  radiative  transfer  in  a  finite  cloud  in 
hypothetical  and  realistic  atmospheres.  Sensitivity  tests  were  car¬ 
ried  out  in  the  infrared  10  urn  window  and  the  50-60  GHz  microwave  re 
gions  and  results  were  compared  with  previous  research.  In  addition 
actual  atmospheric  case  studies  were  selected  from  two  days  and  com¬ 
parisons  made  among  the  calculated  brightness  temperature  values  gen 
eratad  by  the  cylinder  model  and  a  plane-parallel  model,  and  the 
brightness  temperatures  observed  by  the  SSM/T  microwave  sounder. 
Finally,  temperature  retrievals  were  performed  using  the  AFGWC  sta¬ 
tistical  method  and  these  retrieved  profiles  were  compared  to  radio¬ 
sonde  observations. 

In  the  development  of  the  cylinder  model  it  was  found  that  due 
to  the  properties  of  the  modified  Bessel  functions,  higher  order 
spherical  harmonic  approximations  to  the  radiative  transfer  equation 
could  be  generated.  This  may  be  an  advantage  over  the  finite  cubic 
and  rectangular  cloud  models,  since  it  may  be  difficult  to  obtain 
approximations  of  greater  accuracy  than  first  order. 

Comparisons  among  cubic,  plane-parallel  and  cylinder  cloud 
models  in  the  10  cm  infrared  window  were  made.  Emergent  intensities 
at  the  cloud  top  and  bottom  were  presented  and  it  was  noted  that  the 
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finite  cloud  models  produced  zenith  and  nadir  solutions  that  were 
larger  and  smaller,  respectively,  then  the  plane-parallel  model. 

Next,  a  weighted  average  intensity  algorithm  was  constructed  to 
approximate  an  observed  value  at  a  far  distance.  Key  points  gathered 
from  this  analysis  included  the  effects  of  the  order  of  approximation 
employed  in  the  various  models  and  the  manner  in  which  multiple  scat¬ 
tering  and  emission  within  a  finite  cloud  manifest  themselves  in  in¬ 
creasing  amounts  of  energy  emerging  from  the  cloud  sides  as  the  cloud 
optical  depth  is  increased.  Finally,  an  analysis  of  cloud  emissivity 
indicated  that  optically  thick  finite  clouds  in  the  10  um  region  may 
not  be  considered  as  black  clouds  because  their  emissivity  values 
approach  a  value  less  then  1.  This  indicated  that  cloud  top  tempera¬ 
tures  of  finite  clouds  may  be  inferred  to  be  too  large  in  the  10  am 
window  region. 

•  Calculations  were  also  performed  in  a  realistic  atmosphere  by 
placing  the  cylinder  cloud  model  in  a  plane-parallel  absorbing  and 
emitting  atmosphere.  Comparisons  were  made  with  the  plane-parallel 
model  in  the  DMSP  SSM/T  microwave  frequencies.  Both  cloud  models 
behaved  in  a  similar  fashion  although  there  were  some  differences  in 
response  noted,  especially  at  smaller  optical  depths.  A  weighted 
average  intensity  study  in  the  microwave  region  was  conducted.  As 
expected  from  the  results  seen  in  the  infrared  study,  the  weighted 
average  intensity  leaving  the  cloud  at  observation  angles  greater 
than  60°  increased  strongly  as  the  cloud  optical  depth  increased. 

This  was  again  attributed  to  the  increased  effects  of  multiple 
scattering  and  emission  within  the  cloud.  Also,  multiple  scattering 
effects  in  the  SSM/T  frequencies  were  shown  to  be  significant. 
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Continuing,  an  analysis  of  the  cloud  solution  with  a  zero  incident 
intensity  distribution  was  carried  out.  It  was  found  that  the  cylin¬ 
drical  cloud,in  this  case,  emits  isotropically.  Finally,  cloud  emis- 
sivity  at  50.50  GHz  was  presented  versus  cloud  optical  depth.  Al¬ 
though  the  same  trends  were  exhibited  as  were  seen  in  the  infrared 
study  the  cloud  emissivity  was  lower  in  value  for  all  approximations 
in  the  microwave  region. 

Case  studies  were  performed  utilizing  DMSP  SSM/T  sounder  data 
in  conjunction  with  standard  meteorological  data.  Cases  were  divided 
into  clear  and  cloudy/precipitating  studies.  In  the  clear  atmosphere 
cases,  two  methods  were  presented  for  estimating  the  surface  emis¬ 
sivity.  A  method  utilizing  the  window  channel  was  employed  to  esti¬ 
mate  the  emissivity,  and  clear  column  calculations  were  performed 
and  compared  to  observed  values  of  brightness  temperature.  With  the 
exception  of  53.20  GHz,  calculated  values  of  brightness  temperature 
for  clear  cases  were  shown  to  be  larger  than  observed  values.  This 
may  be  attributed  to  a  systematic  error  in  the  source  function  as 
it  is  applied  to  the  remote  sensing  problem.  The  53.20  GHz  results 
are  considerably  influenced  by  the  Earth's  surface.  It  was  conjec¬ 
tured  that  the  surface  temperature  was  underestimated  since  tempera¬ 
tures  reported  at  standard  meteorological  instrument  shelter  heights 
are  not  representati ve  of  the  radiating  temperature  of  the  Earth's 
surface.  This  surface  effect  error  may  balance  the  source  function 
error  in  surface-sensitive  channels. 

Cloudy/precipitating  cases  were  compared  with  cylinder  and 
plane-parallel  model  calculations.  Calculations  in  cloudy  conditions 
were  used  to  generate  scatter  diagrams  that  were  similar  to  the  clear 
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cases  and  it  was  deduced  that  non-precipitating  clouds  have  very 
little  effect  on  microwave  radiative  transfer.  Calculations  in  pre¬ 
cipitating  conditions  were  used  to  generate  scatter  diagrams  that 
demonstrated  the  advantages  of  the  cylinder  model  over  the  plane- 
parallel  model.  Results  for  both  models  were  reasonable  for  cases 
with  rainfall  rates  estimated  at  1  mm/hr.  For  the  two  cases  with  a 
rainfall  rate  estimated  at  5  mm/hr,  the  plane-parallel  model  failed 
while  the  cylinder  model  was  quite  reasonable.  These  differences 
were  attributed  to  the  model's  geometries.  It  was  interesting  to 
note  that  the  greatest  failure  of  the  plane-parallel  model  happened 
when  the  precipitating  cloud  system  modeled  was  most  cylinder  like. 

Finally,  an  exercise  in  temperature  retrieval  was  performed 
using  the  AFGWC  statistical  method.  Here  it  was  demonstrated  that 
precipitating  conditions  cause  a  consistent  deviation  of  the  retrieved 
temperature  profile  from  the  observed  radiosonde  temperature  profile 
in  the  lower  boundary  layer  -  as  much  as  15°  K  in  5  mm/hr  rainfall 
rates  covering  only  30%  of  the  field  of  view. 


APPENDIX 


To  carry  out  the  operation  indicated  in  Eqs.  (5-18)  and  (5-58), 
the  expansions  for  the  intensity,  Eq.  (5-11),  and  for  the  phase 
function,  Eq.  (5-10),  are  substituted  into  Eq.  (5-17).  In  this  in¬ 
stance,  however,  we  use  the  two-dimensional  representation  for  in¬ 
tensity,  I(r,^).  Multiplying  by  the  complex  conjugate  of  the  spheri¬ 
cal  harmonic  and  integrating,  Eq.  (5-17)  may  be  represented  as 
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The  integrations  on  the  right  hand  side  of  Eq.  (A-l)  may  be  done 
straightforwardly.  To  accomplish  the  left  hand  side,  the  recursion 
formula,  Eqs.  (5-15)  and  (5-16),  may  be  added  to  produce  the  term 
sinecoscuY^^) ,  and  subtracted  to  produce  sinesincuY^^) .  The  details 
of  the  sin9sin<|>Y™(£)  are  presented  here.  First,  the  constant  of  in¬ 
tegration  from  the  right  hand  side  is  multiplied  through  the  equation, 
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which  may  be  written  as 


VITA 


Name 

Bi rthdate 
Bi rthplace 
High  School 

Un  i  ve  rs  i  ti  es 


Degrees 


Professional  and  Honorary 
Societies/Memberships 

Publi  cati  ons 


Grant  Clifford  Aufderhaar 

August  20,  1948 

Celina,  Ohio 

Memorial  High  School 
St.  Mary's,  Ohio 

United  States  Air  Force  Academy 

Colorado 

1966-1970 

Texas  A  &  M  University 
College  Station,  Texas 
1971-1972 

University  of  Utah 
Salt  Lake  City,  Utah 
1977-1980 

8.  S.,  1970 

United  States  Air  Force  Academy 
Colorado 

M.  S.,  1972 

Texas  A  &  M  University 
College  Station,  Texas 

American  Meteorological  Society 
CHI  Epsilon  PI 
PHI  Kappa  PHI 

"Remote  Sounding  of  Cloud  Composi¬ 
tions  from  NQAA  IV  and  Nimbus  VI 
Infrared  Sounders,"  1977  with  K.  N. 
Liou,  R.  G.  Feddes  and  T.  L. 

Stoffel .  AFGL-TR-77-0252 . 


1/2  (o,-6*2)’/25|-,f:'  (r.t) 


d<Pr  ct+1  v  ,(>r 


Part  B  is 


V  m^T)  (cH-S-l)172  (af6)1/2^-lS-|  (r,^} 


Part  C  is 


-  +  ,  1  1 

Ts-F7C%r- 


(=.-S-l)V2  (c.-e> 1/2  rf  l£>  (r ,<p  ) 


34>r  a-1  '  ,cpr' 


Part  0  is 


JL_  yS+1  /  . 

9<J>  a+1  '  ,:V 


The  second  term  in  Eq.  (A-4)  also  has  four  portions  which 
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Combining  the  first  and  second  terms  in  Eq.  (A-4),  and  the  sinecos-j 

terms,  and  defining  A®,  8®,  C8,  D8  as  before,  we  arrive  at  Eq.  (5-59) 
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Note  that  if  the  intensity,  I(r,2)  is  reduced  to  the  one-dimensional 


case,  I(r,n),  then  the  term  is  zero  and  Eq.  (5-59)  reduces  to 
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Eq.  (5-19). 
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